NASA TECHNICAL 
REPORT- 





Ft 

COPY 



\ 

FOURIER TRANSFORM REPRESENTATION 
OF AN IDEAL LENS IN COHERENT 
OPTICAL SYSTEMS 

by Gerald J. Grebowsky 

Goddard Space Flight Center 
Greenbelt, Md. 

NATIONAL AERONAUTICS AND SPACE ADMINISTRATION * WASHINGTON, D. C. • 



IR R-319 


JANUARY 1970 



NASA TR R-319 


FOURIER TRANSFORM REPRESENTATION OF AN 
IDEAL LENS IN COHERENT OPTICAL SYSTEMS 


By Gerald J. Grebowsky 

Goddard Space Flight Center 
Greenbelt, Md. 


NATIONAL AERONAUTICS AND SPACE ADMINISTRATION 


For sale by the Clearinghouse for Federal Scientific and Technical Information 
Springfield, Virginia 22151 - Price $3.00 



ABSTRACT 


This document presents a mathematical analysis of the 
approximations required to obtain the Fourier transform 
representation of an ideal lens. An attempt is made 
throughout the paper to demonstrate the physical signifi- 
cance of the approximations, and the variations from ideal 
results, produced by neglected terms in the mathematical 
formulation. The approximations involved are considered 
in terms of the output signals in optical spectrum analyzer, 
optical imaging, and optical correlator systems. 
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FOURIER TRANSFORM REPRESENTATION OF AN 
IDEAL LENS IN COHERENT OPTICAL SYSTEMS 


t>y 

Gerald J. Grebowsky 
Goddard Space Flight Center 


INTRODUCTION 

In recent years there has been a growing interest in the application of optical imaging tech- 
niques for the purpose of processing data signals. These efforts are largely based on the interpre- 
tation of optical imaging systems as spatial filters (Reference 1). By introducing Fourier transform 
methods the relation between an object and its image has the same form as the relation between 
the input and output signals of an electrical system. Comparing the transform of the object and 
image, the imaging process (unity magnification is assumed here) can be described by 

!K' u y) = T K* ■“») • 


or 


Image spectrum - Transfer function x Obj ect spectrum . 


This expression has the same form as that of an electrical network except that the optical spec- 
trums are two-dimensional. Since optical objects and images are two-dimensional, a Fourier 
transform must be taken with respect to two spatial coordinates instead of the single time coordi- 
nate which appears in electrical signals. 

It is the transfer function relation given above that leads to the spatial filter interpretation of 
optical imaging systems. The optical transfer function, T^ x , wj, is a characteristic of the optical 
elements in a system. An ideal imaging system would have a transfer function that is constant 
over the object frequency range. In such an ideal system the image would be an exact replica of 
the object. The corresponding electrical system would have a flat frequency response over the 
bandwidth of the input signal. 

It may appear at first thought that the transfer-function notation is nothing more than an 
arbitrary selection of notation. However, in optical imaging systems, an object represented by a 
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sinusoidal spatial variation of light amplitude is imaged as such— even in the presence of aberra- 
tions (Reference 2). Aberration effects appear as reduced contrast and lateral shift of the sinus- 
oidal test image. Thus, using sinusoidal test gratings, it is theoretically possible to experimentally 
determine the transfer function for a given optical system, hi general, the optical transfer func- 
tion can have complex values. The magnitude is related to the reduction in contrast, and the phase 
is related to the lateral shift of the image. In an actual system, the transfer function will not have 
the constant amplitude and phase of the ideal imaging system described above. 

Since the optical transfer function is determined by comparing the light amplitude variations 
of the output image to that of the original object input, the introduction of any additional element 
into the optical system to vary the amplitude and/or phase transmission properties will change the 
optical transfer function. To readily utilize an optical system as a spatial filter in a fairly direct 
manner, it is necessary to know what amplitude and phase variations should be inserted and where 
they should be inserted. Otherwise, obtaining a particular transfer function for a spatial filter 
application would be a trial and error proposition. This implementation problem is solved for 
many cases of practical importance by the optical Fourier transform representation discussed in 
this report. 

Within certain limitations, the light amplitude distribution in the back focal plane of a lense is 
proportional to the two-dimensional Fourier transform of the light amplitude distribution of a two- 
dimensional object inserted on the front side of the lens. Within the range of validity for this opti- 
cal Fourier transform representation, the transfer function is varied by a multiplicative factor 
represented by the amplitude and phase transmission properties of an element inserted into the 
back focal plane of the lens. For example, to set the transfer function equal to zero for a particu- 
lar frequency component, the light passing through the corresponding point in the back focal plane 
of the lens is simply blocked. 

The mathematical development of the optical Fourier transform representation presented in 
this report is intended to clarify the limitations and interpretation of the Fourier transform opera- 
tion of a lens. The derivation is based on the Rayleigh- Sommerfeld diffraction formula and on 
optical paths defined by geometrical ray tracing. As each limitation is introduced, an attempt is 
made to describe the effects on the accuracy of the optical Fourier transform representation. Such 
detailed consideration has been found lacking in available treatments of the derivation (Reference 3) 
and is the main purpose for the development presented in this report. 


FOCAL PROPERTIES OF A LENS 

To derive the formula for a focussed diffraction pattern, the properties of an ideal lens are 
defined. Our discussion is restricted to the case of an ideal lens; the effects of lens aberrations 
and diffraction at the edge of the lens are ignored. These effects are assumed to be taken into ac- 
count by modifying the end result or by restricting the range of variables to a region in which the 
ideal assumptions are valid within experimental accuracies. 
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The definition of an ideal lens is based on the geometrical focussing properties shown in Fig- 
ure 1. The properties assumed can be stated as follows: 

1. The lens can be represented by a plane, l, perpendicular to the optical axis, with all re- 
fraction taking place at this plane. This is a thin lens approximation that neglects the 
thickness of the lens. 

2. The rays passing through the point o (intersection of the optical axis and the lens plane l) 
are called principal rays; they are not deviated. 

3. All incident rays parallel to a principal ray will be focussed at the point where the prin- 

cipal ray intersects the back focal plane F'. That is, the light reaching a point in the back 
focal plane, F', at a distance, p - f tan e (where f is the focal length— the distance between 
the planes L and F') from the optical axis, is contributed by a principal ray making an 
angle 9 with the optical axis plus all rays initially parallel to this principal ray. ‘ 

4. If a plane p is constructed perpendicular to a bundle of parallel incident rays, the optical 
path length will be the same for any of the parallel rays moving from p to the common 
point of focus in F'. 

For any set of parallel rays, Figure 1 represents the projection of the parallel rays onto the 
plane through the optical axis and the principal ray as shown in Figure 2. In Figure 1 the distances 
p and p x are measured from the optical axis in the plane of the figure. As shown in Figure 2, these 
distances represent different quantities in the planes F and F' respectively. The distance p in the 
back focal plane F' is the distance from the optical axis to a point (x, y). In the plane F, p is not 
the distance from the optical axis to a point ( x !> y x ) ; P x is the projection of this distance onto the 
axis defined by the intersection of the plane Fwith the plane through the optical axis and a point (x,y) 
in the plane F'. Since the orientation of the p 1 axis will depend on the position (x, y ) in F', p will be a 
function of x, y, x , and y whereas p depends only on x and y. The difference in meaning of p and p 1 
also appears when the algebraic sign is considered. In Figures 1 and 2, p is the distance from the 
optical axis to the point of focus and is a positive quantity regardless of where the point is located. 
On the other hand, p x is a coordinate of the intersection of a ray with the F plane and we will use 


F P L F' 



Figure 1— Ideal focusing of parallel rays. 


Figure 2— Difference between p and P 1 


3 




(a) F' PLANE 

y i 



the sign convention that p 1 is positive above the 
optical axis and negative below it when drawn as in 
Figure 1 (rays sloping down to the right). In refer- 
ence to a point (x, y) in F', the positive p 1 axis will 
lie in the quadrant opposite to that of the point (x, y). 

Figure 3 shows the geometry of the various lines 
in planes F and F' and demonstrates the sign conven- 
tion for p v From the geometry of Figure 3(a), the 
following relations are found: 


, _ X y 

cos <p - — > sin 4> - — > 

P = ( X 2+y2)l/2. (1) 

From the geometry of Figure 3(b), p can be deter- 
mined for any point ( x x , y x ) as follows: 


<t> ~ a + / 3 , 


cos a 



r 



r = (Xj2 + y2)l/2 _ 


Pi - - r cos/3 - - r cos (0“a) 

= - r {cos cos a + sin<£ sin a} » 


or 


Figure 3— Geometry for p and p 1 . P i x i cos0 y 1 sin^ . 

Substituting the relations for cos <j> and sin^ obtained from Figure 3(a), the expression for p can 
be rewritten as 


- xx - yy - XX ' yy 

P \ ~ p “ ( x 2 + y 2 ^ 1/2 ' (2) 

The dependence of p on x and y and of p l on x, y, x , and y is explicit in Equations 1 and 2 respec- 
tively. These results agree with the discussion in the last paragraph. 

Figure 4 is now used to determine an expression for the optical path length from a point (x x , y x ) 
in F to a point (x, y) in F' . Figure 4 shows a principal ray AA' and a representative parallel ray 
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(dotted). The plane P is perpendicular to the 
incident rays. Since the thickness of the lens 
is neglected, the optical path length r for the 
principal ray AA' will be given by the geomet- 
rical length: 


r 


+ £ 


(3) 


From the geometry of Figure 4, l x is given by 
the expression 



Figure 4— Geometry for optical path length. 


- (f 2 +,£>2)1/2. 


By similar triangles, l is given by 


2 

d 


f 


df 


(4) 


(5) 


For the principal ray, AA', the optical path length from the plane P to the point A' in the back focal 
plane F' is + ^ 2 . Using Equations 4 and 5 this length is given as 




2 


df 


if + df 


f 2 + p 2 + df _ 

(f 2 + p2)V2 ‘ 


By the fourth assumption for an ideal lens, the optical path length from the plane p to the focus 
point A' is the same for all the rays parallel to AA' (note that l +^ 2 does not depend on p x ). There- 
fore the expression found for + ^ 2 holds for every parallel ray and the expression for the optical 

path length r can be written 


r 


^ + 't n 


f 2 + df , 

(f 2 +p 2 jl /2 S 


( 6 ) 


The term 1 2 remains to be determined. By comparison for the principal ray and the representative 
ray in Figure 4 it should be obvious that this term will not be the same for all parallel rays. From 
the right triangle with 1 2 as a leg, 1 3 can be expressed as 

= p x sin 0 . 
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However, sin e - p/[ f 2 + p 2 ) V2 , so the expression for 1 3 can be rewritten as 


PP 


(f 2 + P 2 ) 


2\V2 


Substituting for l 3 in the expression for the optical path length r ? then 


( 7 ) 


f 2 + df + p 2 + pp^ 

T ~ (f2 ~p 2)1/2 (8) 

Previously derived expressions for p and p are given by Equations 1 and 2. Substituting for p 
and pj in Equation 8 we obtain 

f 2 + x 2 + y 2 + df - xx 1 - yy 1 

r ~ (f 2 + X 2 + y 2 ) 1/2 ( 9 ) 

This expression gives the optical path length from any point (x x , y i ) in a plane F (a distance d in 
front of the lens) to a point (x, y) in the back focal plane F'. To facilitate further discussion, this 
expression can be written 


r - R(x, y) - aXj - fiy t , 


(10) 


where 


f 2 + df + x 2 + y 2 
R(X ’ y) (f 2 +x 2 +y 2 ) 1/2 

X 

a “ (f 2 + x 2 +y 2 )!/ 2 ’ 



Figure 5— Focussing of a cone of light. 


and 

B = 7 2 

(f 2 +X 2 +y 2 )l/2 

To summarize our results at this point, 
consider a point A in the plane f as shown in 
Figure 5. Assuming that light is radiated in 
all directions from the point A, consider the 
portion of light propagated in directions at an 
angle 6 with respect to the normal, n, to the 
plane F. The light rays representing these 
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directions will form the surface of a cone of half angle 6 as shown in Figure 5. For each of these 
rays a parallel principal ray can be drawn and each principal ray will make an angle 6 with the 
optical axis. Thus for each ray at an angle 0 with respect to the normal to F at the point A, the 
results derived above can be applied. That is, each ray is focussed to a point on the ring of radius 
P = f tan# where the corresponding principal ray intersects the back focal plane F'. The optical 
path length from A to each point on the ring is given by Equation 10. Since this holds for any point 
A, the general focal properties of the ideal lens can be stated: 

1. The light radiated from all points on F in directions at an angle 0 with respect to the nor- 
mal to F is focussed into a ring of radius p - f tan 6 in the back focal plane F' . 

2. The optical path length r from any point (x i , y i ) on the plane F to any point (x, y) in the 
back focal plane F' is given by Equation 10. 

In these statements the plane Fisa plane perpendicular to the optical axis at a distance d in front 
of the lens. 

Anticipating the derivations in the next section, the relation between p and 6 , specified by the 
first focal property above, can also be expressed in terms of cos 0 as 

„ _ f _ f 

(f 1 2 + p 2 ) 1/2 (f 2 + x 2 + y 2 ) 1/2 (11) 

This expression can be derived from the geometry of the figures in this section or from the equa- 
tion, p - f tan 9 , in terms of tan 0 by applying trigonometric identities. 


FOCUSED DIFFRACTION PATTERN 

Since only light distributions on plane surfaces are being considered, the Rayleigh- Sommerfe Id 
diffraction formula can be used (see Appendix A). In rectangular coordinates this diffraction for- 
mula has the form: 


A(x, y, z) 



cos 6 dx x dy x . 


( 12 ) 


This formula gives the complex light amplitude A(x, y, z) at any point in space (z > 0) caused by a 
monochromatic coherent light distribution A' (x x , y x ) given for every point (x x , y x ) in a plane f. 
Referring to Figure 6, the terms in the diffraction formula are defined as 

1. A' ( x x , y x ) is the complex amplitude of monochromatic light given for all points ( x x , y x ) in 
a plane F located at z = 0. 

2. A(x, y, z) is the complex amplitude of light produced by A(x x , y x )at a point (x, y, z) in 
space (z> 0). 
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Figure 6— Relation between points (x, y, z) and 
( x ]7 y 1 ) in diffraction formula. 


3. r is the distance from a point y x ) 
in plane F to the point (x, y, z). 

4. e is the angle between r and n, where 

r is directed from (x x , y x ) to (x, y, z) 
and n is the normal to the plane f at 
( x , y i ) in the direction of the positive 
z axis. The term cos 6 is usually re- 
ferred to as the obliquity factor. 

5. k = 2tt/\ , where k is the wavelength of 
the monochromatic light. 

In general each of the three vector components 
of the electromagnetic field representing the 
light distribution must be determined by the 
diffraction formula. This discussion considers 
the scalar light amplitude distribution which 
requires only one equation (Reference 4). hi 


practice this is permissible if polarization affects can be neglected. Thus the light amplitude can 
be defined such that the square of its absolute magnitude gives the intensity, which is a measurable 
quantity. 

The diffraction formula can be immediately simplified by considering the relative magnitudes 
of the terms inside the brackets: 



For wavelengths k as long as 100 microns (far infrared) the first term is relatively large (600) 
while for r larger than 1 cm the second term is less than one. For visible light \ is much less 
than 100 microns (0.4 to 0.7 microns), and k is of the order of 10 s . Since this discussion (as in 
most cases in optics) will deal only with r greater than one centimeter, the 1/r term is negligible 
and can be dropped without any appreciable effect on accuracy. The diffraction formula, Equation 
12, can therefore be written 


A(x, y > z) 



co s 0 dx l dy 


(13) 


where the constant factor, ik, has been taken outside the integral. 

The- obliquity factor, cos 9, is a weighting factor that accounts for the difference in the amount 
of light radiated in different directions. Since cos 9 has a maximum value of one when 9 equals 
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zero, this factor has a maximum value of one 
for light contributions propagated normal to 
the signal plane and drops off as the angle with 
respect to the surface normal increases. Re- 
ferring to Figure 6, if the light from a point 
( x r y i ) contributing to the light at the point 
(x, y, z) is assumed to travel the straight line 
r, this line is a light ray at an angle 9 to the 
normal n. In the previous section, it was 
shown that through the focal property of an 
ideal lens, this angle is a constant for all light 
contributing to a point (x, y) in the back focal 
plane and that x, y, and 9 are related by the 
expression 


C ° S(? = ( f2+x 2 f +y 2)l/2 • (11) 

In other words, an ideal lens focuses light of 
constant obliquity factor into a ring of radius 
( x 2 + y 2 )i /2 specified for a given 9 by the above 
expression. 


F PLANE F' PLANE 




The significance of this focal effect can be _ _ . ...... . . 

Figure 7— Comparison of diffraction with 

seen by comparing the two diagrams in Fig- and w!t h 0 ut focussing, 

ure 7. In Figure 7(a), the points A and b are 

sample points in the F plane, and the points c and D are sample points in a parallel plane at z = d + f . 
Considering the point c, it is noted that the paths AC and BC have obliquity factors of cos 6 and 
cos d 2 respectively. From this example it is obvious that for a point such as c (or d) the obliquity 
factor will depend on the location of the contributing point (x x> y x ) in F. Likewise, considering the 
point A in f, it is noted that the paths Ac and AD have obliquity factors of cos e and e 3 respectively. 
This indicates that the obliquity factor also depends on the location of the point (x, y, z). Since de- 
termining the obliquity factor is included in the derivation of the diffraction formula, it is given 
here without proof for the case z = d + f as shown in Figure 7 (a) : 


cos 0 


d + f 

[(x- x x ) 2 + (y- y i) 2 + (d+ f) 2 


(14) 


This expression includes the coordinates of both the point ( , y x ) , in the source plane, and the 
point (x, y, z), at which the diffracted light amplitude is to be found. The expression for .the general 
case will have a z in place of the (d + f) used for the special case of Figure 7(a). In the diffraction 
formula, the obliquity factor appears under the integral sign since x x and y x are the variables of 
integration and appear in the obliquity factor as given by Equation 14. 
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Let us now consider the case of focussed diffraction as illustrated in Figure 7(b). In Fig- 
ure 7(b), only the rays AC and BD of Figure 7(a) are considered, and as indicated by the angle 0 , 

AC and BD are parallel rays. The dotted portion of these rays indicates the path of light followed 
in Figure 7(a). Because of refraction by the lens, these paths are changed to those focussed to the 
point E. Now, considering a point such as E, then the obliquity factor cos 6 x is the same for points 
A and B and therefore independent of the coordinates (x^ y x ) of the point in F. Considering any 
other point G in F' , it is recalled that to contribute to a point G a ray must be parallel to the principal 
ray OG. Rays parallel to OG will have an obliquity factor cos 6 different from cos 6 for the point 
E. Thus the obliquity factor does depend on the location of the point (x, y) in the back focal plane. 

The obliquity factor for the case of a focussed diffraction pattern is given by Equation 11 as 
cos Q = £/( f 2 + x 2 + y 2 ) 1/2 and does not depend on the coordinates x and y . 

Since the obliquity factor for the focussed diffraction pattern is independent of the integration 
variables x t and y , this factor can be taken outside the integral and Equation 13 can be written as 


AO, y) 


\(f 2 




r 


dxj dy j - 


(15) 


where A(x, y) now represents the complex light amplitude at a point (x, y) in the back focal plane 
of a lens. 

To complete the discussion, now consider the term r which was defined as the distance from 
the contributing point to the point of interest. In Figure 6 this distance is measured along the 
straight line from ( x x , y x ) to (x, y, z). In Figure 7(b), the light traveling from A to e does not fol- 
low a straight line because of refraction at the lens plane L . Assuming that the effects of the length 
of the refracted path are the same as traveling an equivalent distance in a straight line, then the r 
in the diffraction formula can be interpreted as the optical path length determined in the previous 
section. Thus e ikr represents the change in phase over an optical length r, and 1/r is an attenua- 
tion factor that decreases the amplitude contribution as the optical path length increases. 

Substituting the optical path length expression for r, as given by Equations 9 and 10, Equa- 
tion 15 becomes 


A(x, y) 


q 


_i jp /*/• A 

f 2 + x 2 + y 2 ) 1/2 f“£ 


A' (x 2 , yj e 


ikf<x,y)-ax j-/3y J 


F 2 + x 2 + y 2 + df - xx 1 - yy 1 


- dXj dy j 


(f 2 + x 2 + y 2 ) 


1/2 


(16) 


The terms (f 2 +x 2 + y 2 )i /2 cancel and e ikR < x ' y > can be taken outside the integral to give 


A(x, y) 


-i± e ikR(x,y ) JJ 


“i27T [ax j/\+/3y 


A< K- y t ) e 

f 2 + x 2 + y 2 + df - xXj - yy^ 


dXj dy x . 


(17) 
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Now introducing the new variables p and q defined as 


a_ _ x 

A " \(f2 + X 2 +y2)l/2 


and 


_ £_ _ y , 

^ A A.(f 2 +x 2 +y 2 ) 1/ 2 

and factor f 2 +df from the denominator of the integral of Equation 17, then 


(18a) 


(18b) 


A(x, y) 


\(f + d) 


; ikR(x,y) Jj" 


>(*- *1) * y(y~yO 


dXj dy t 


1 + 


f(f + d) 


( 19 ). 


FOURIER TRANSFORM APPROXIMATION 


The algebraic identity, 


can be used to obtain 


l + 


l - 


1 + I 


1 

x(x- x x ) + y(y- y x ) 

F(TTd) 


1 

HI + d ) 

x(x-x 1 ) + y(y-y 1 ) 


Introducing this identity, Equation 19 can be written as 


A ( x > y) - \(f + d) 


*ikR(x,y) 


JJ 


l- 


i + 


f(f + d) 


A' (x x , y x ) 


£ ( x_ x i) + y(y~yiU 

Equation 20 can be rewritten with an integral for each term in the bracket so that 


( 20 ) 


A(x, y) - A.(f+d) 


a ikR(x,y ) 


JJ A' ( x ,. y x ) 


‘ i27r ( PXl+qy O dXj dy 1 


T \( f .+ d) 


i ikR(x t y) 


•r A 

JJ i 


A' (x , y x ) 


-^(pXj + qy ^ 


F(f + d) dx i d yi • (2i) 


x(x-x x ) + y(y-y x ) 
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By restricting the maximum values (aperture limits) of x, y , x ] , and the second integral of 
Equation 21 can be made negligible compared to the first since the denominator of the integrand can 
be made large. This approximation will be discussed in more detail later; here it is simply as- 
sumed that it is possible to neglect the second integral. 

The diffraction formula can then be written approximately as 


-i e ikR(X,y) 

A ( x - y) = \(f + d) f (p> <0 ’ 

where F(p, q) is the two-dimensional Fourier transform of A' (x , y } given as 


F(p, q) 



(PXi+qy,) 


dXj dy j 


(22) 


To measure this light distribution, it is the intensity which is the square of the magnitude of the 
complex amplitude A(x, y), which is measured. This intensity is given as 


I(x, y) 


A(x, y) A* (x, y) 


1 

[A( f + d)] 2 


F(p, q) F* (p, q) . 


(23) 


or 


i( x > y) 


lF(p, q) 1 2 _ 
[A(f + d)] 2 


Thus the intensity in the back focal plane of a lens (within the limits to be determined for the ap- 
proximation made) is given by the square of the magnitude of the Fourier transform of the light 
amplitude in the plane F. 

If the aperture restrictions in the back focal plane limit the maximum values of x and y so that 
the phase variations caused by the exponential term in front of the integral in Equation 22 can be 
considered constant, then 


A(x, y) - KF(p, q) , 


(24) 


where K is a complex constant given by 


-ie ikR(x,y) 

\(f + d) ‘ 
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Thus, within the range of x and y for which e ikR ( x * y ) can be assumed constant (i.e. negligible phase 
variation), the amplitude distribution in the back focal plane is proportional to the Fourier trans- 
form of the light amplitude distribution in the plane F. This relationship requires tighter restric- 
tions on x and y than our previous approximation. In terms of spectrum analysis in the back focal 
plane, this relation is not important since only intensity can be measured. However, in cascaded 
lens systems, the Fourier transform relation between amplitudes allows each pair of lenses to be 
accounted for by a double Fourier transform operation. The advantages of such a relation will be 
demonstrated in a later section. 

FOURIER COMPONENTS (p, q) AND FOCAL PLANE COORDINATES (x, yj 

In the previous discussions, the amplitude distribution A(x, y) has been expressed in terms of 
the Fourier transform of A' ( x , y ) (refer to Equation 22). However, the Fourier transform co- 
ordinates are p and q and are defined by Equations 18 as 


X 

p ~ \(f* +X 2 + y 2 ) 1/2 ’ (18a) 

_ y 

q \(f 2 + X 2 +y2)l/2 ‘ (18b) 

Substituting for p and q in Equation 22, an expression for A(x, y) can be obtained in terms of x and 
y; however, this result is somewhat complicated by the fact that p and q are each dependent on 
both x and y. When the above expressions are substituted for p and q, the optical Fourier transfori 
is given by 


F(p» q) - F 


x , y \ 

A. (f 2 + x 2 + y 2 ) 172 \ (f 2 + x 2 + y 2 ) 172 / 


JJ A ' ( x i’ yJ e 1 


2 77 [xx 2 +x 2 +y 2 ) 1 7 2 + y y jA(f 2 + x 2 + y 2 ) 1/2 ] 


dx x dy x 


(25) 


A more desirable relationship would exist if p were directly proportional to x and if q were directb 
proportional to y. Then the light amplitude at a particular value of x would be related to a particuli 
value of p and a similar relation would exist between q and y. As given by Equation 18, p.and x 
(also q and y) are not so simply related since p also depends on y (q also depends on x). In a later 
discussion the importance of a linear relation between the transform coordinates, p, q, and the 
spatial coordinates, x, y, will be shown. 
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To obtain a linear relation between p and x, consider the series expansion of Equation 18a. 


x 

x 2 + y 2 3 /x 2 + y 2 \ 2 

p " \f 

} 2f 2 1 8 \ f 2 / 


Restricting our analysis to an area of the focal plane such that 

x 2 + y 2 

£_ « i , 

2f 2 

.then all but the first term in the brackets can be neglected to obtain the approximation 


P = XT’ (26) 

Similarly, an approximation of q can be obtained 

q = 7T ' (27) 

Thus, within a restricted area of the back focal plane, the Fourier transform expression can be 
written as 


F(p. q) 


= JJ A '( x i’ 


y y e -i2TT(x Xl Af+y yi Af) 


dXj dy j • 


(28) 


The coordinates x and y in the back focal plane are then scaled representations of the frequencies 
p and q respectively. That is, light contributions corresponding to a spatial frequency p in the 
direction appear at the coordinate x = X.fp in the back focal plane. (Similarly, contributions of 
spatial frequency q in the y x direction appear at y = \fq.) 

The actual restriction to be imposed on x and y for the above approximation will depend on 
how accurate a Fourier frequency value is required in a particular application. The error in the 
approximate frequency of Equations 26 and 27 as a fraction of the exact value given by Equation 
18 is 


(xr) 

X ' 

\(f 2 + x 2 + y 2 ) 172 , 



x 2 + y 2 
f 2 



1 . 


(29) 
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Letting r 2 = x 2 + y 2 (r is the radius of a circle in the x, y plane), then r can be expressed in terms 
of multiples of the focal length f such that 


r af . 


( 30 ) 


Substituting r 2 = x 2 +y 2 = a 2 f 2 into Equation 29, 


E f = ( 1 + a 2 ) 1/2 - l . 

The curve of Figure 8 gives the percent error (100 E f ) of the linear approximation of frequency as 
a function of a. For a less than 0.14, the error will be less than 1 percent. Thus the linear 
approximations 


p = XT' (26) 

and 

<* = XF ’ (27) 

are accurate to within 1 percent for values of x 
and y satisfying the restriction 

(x 2 +y 2 ) 1/2 = r ~ 0 . 14 f . 

For accuracies better than 1 percent, smaller 
values of a must be imposed as given by Equa- 
tion 3 1 and Figure 8. 



r 

a = 7 

Figure 8— Percent error in 
linear frequency approximation. 


Since the approximation requires limiting consideration to the area within a circle of radius 
r m ax = 0.14f (for accuracy within 1 percent), the maximum value of x 2 +y 2 is specified by 


(x 2 +y 2 ) n 


0. 02f 2 


Squaring the approximate expressions for p and q (Equations 26 and 27) and adding, 


p 2 + q 2 


x 2 +y 2 
f 2 X 2 


(32) 
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Applying the restrictions on x 2 +y 2 to Equation 32, 


The maximum allowed value of p occurs when q = 0, and the maximum q occurs when p = 0 

0.14 _ . 

Pmax - “X" for q - 0 ( 1 . e. , y - 0) . 

0.14 

q ma x “IT for P = 0(i.e., X - 0) • 

As an example, consider green light of wavelength 5461 x 10' 8 cm. In this case the simplified 
expressions for p and q are accurate to within 1 percent for frequencies in the range 

/ „ < 0.14 _ o. 14 x io 8 

(p 2 + q 2 ) 1/2 = - y - % 2560 cycles/cm . 

On the x axis (y = 0, q = 0), the maximum frequency will be 

^max 2560 cycles/cm . 

On the y axis (x = 0, p = 0), the maximum frequency is likewise 

q max = 2560 cycle/cm . 

In practice, the limitation of available techniques for controlling the input light distribution, 

A ' ( x i , y i ) , restricts maximum spatial frequencies to values below the 2560 cycles/centimeter re- 
striction imposed above. Therefore, the linear approximations of the frequency components p and 
q (Equations 26 and 27) are applicable for practical systems and Equations 22 and 23 can be ex- 
pressed as 


A(x, y) 


KF 1 


fe ’ \f) 


(32) 


I(x, y) 



[X(f + d)] 2 


(33) 


The accuracy of p and q used in Equations 32 and 33, and as defined by Equation 18, is deter- 
mined by Equation31, as shown in Figure 8 for values of a = r/f . Thus for a given focal length, f , 
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the restriction on the maximum value of r determines the accuracy of the linear approximation 
introduced here. Of course, these equations also include the approximation assumed earlier in 
neglecting terms other than the F(p, q) term. The next section considers that approximation and 
whether the restriction x 2 + y 2 ^ 0.02 f 2 is a sufficient restriction to assure the validity of neglecting 
terms other than F(p, q) . 


APPROXIMATION LIMITS FOR THE FOURIER TRANSFORM REPRESENTATION 

Returning to the focussed diffraction formula given by Equation 20, 


A(x, y) 


_ ie ikR(*.y) ff 

Mf + d)JJ 1 


1 + 


1 

f(f +d) 

x(x-Xi) + y(y - y t ) 


A ' ( x i> *i) e 1 


i2?T(px j + qy j ) 


dx 


i d yi 


(20) 


consider the limitations required to obtain the Fourier transform approximation given by Equa- 
tion 22. To obtain the form of a Fourier transform of A' (x 2 , y 2 ) , the bracketed term must be ap- 
proximated by a constant. This term can be assumed equal to one, if the range of x, y, x 1 , and y 2 
is restricted to satisfy the inequality 


1 + 


1 

f(f + d) 

X 2 + y 2 - (xx 2 +yy 2 ) 


« 1 . 


Referring back to Equation 21, this approximation corresponds to making the second integral neg- 
ligible compared to the first integral, which has the form of a Fourier transform. The complete 
term inside the brackets of Equation 20 is effectively a weighting factor which varies the contri- 
bution from each point (x x , y ) to the point (x, y). This factor represents the effect of the obliquity 
factor and path length attenuation. It is usually assumed that these effects are negligible and that 
the inequality condition is satisfied. In the following analysis, a more detailed quantitative dis- 
cussion of this approximation is presented. 

Neglecting the variable term, when it satisfies the inequality condition given above, is an ap- 
proximation of the light amplitude contribution from each point ( x 2 , y 2 ) to a point (x, y). That is, 
the contribution dA(x, y) at a point (x, y) from an infinitesimal region dx 2 dy 2 about the point ( x 1 , y ) 
is given exactly by 


dA(x, y) - K 


1- 


1 + 


f(f + d) 


: ( x_ x i) + y(y _ yiU 


A' (x x , y 2 ) e 


•i 2^(px 1 + qy 


dx x dy x 


and applying the approximation of neglecting the variable term inside the bracket, then 


dA(x, y) 


KA' (x 1>yi ) e 


^^(pXj + qy j) 


dXj dy x . 


(34) 


17 



The total light amplitude A(x, y) at a point (x, y) is obtained by integrating over the range of x t 
and y . The integration of Equation 34 will yield an approximation for the total light amplitude 
A(x, y) at least as accurate as the worse case of Equation 34. That is, the greatest possible error 
would be given by the maximum value of the neglected term. 

Thus, to determine the limitations to be imposed, consider the maximum error introduced by 
neglecting the variable term to obtain Equation 34. Denoting the error by the fraction, E A , given by 
the ratio of the neglected term to the exact factor within the brackets of Equation 20, then 


x 2 +y 2 - (xxj + yyj 

E a = f(f + d) 


(35) 


To simplify the discussion, x, y, Xj , and y x are expressed in terms of polar coordinates r, <f>, 
r x , and 4 > x . The relations between these coordinates are given by 


x 2 + y 2 , 

n 

Cl 

U 

x t 2 + y x 2 

r cos 0 , 

x i = 

r x cos^> 1 

r sin0 , 

= 

r x sin0 x 


Substituting in Equation 35 


r 2 - ri-j ( cos cos + sin sin 
f(f + d) 


(36) 



Figure 9— Diagram of conditions for maximum E A . 


Using the identity cos (0-^) ~ cos^cos^ 

+ sind> sin <j> 1 , Equation 36 can be rewritten as 


r 2 - rr j cos (<£ *d>j ) 
E a = f(f + d) ' 


(37) 


The cosine term can take values between 
plus and minus one. Since the maximum error 
is of interest, consider the case for cos (d>- 
= - 1. Equation 37 can be rewritten for cosine 
equal to minus one as 


r 2 + rrj 

e a ~ f (f + d) ’ (38) 

Figure 9 shows the relative positions of points 
(x t , y x ) and (x, y) for the case cos = - 1. 
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As shown by the figure, the maximum error defined by Equation 38 applies to the light contribu- 
tions from points (x , y x ) located on the line of intersection 0 1 T between the planes F and Q. The 
plane Q is a plane containing the optical axis and the point (x, y) in the back focal plane F' . The 
points (x lf y i ) are further restricted to the portion of the line of intersection of F and q that is on 
the side of the optical axis opposite from the point (x, y). From the geometry of the figure it is 
clear that the angle <t> is equal to tt + <j> l . Thus <t> ~ <£• is equal to tt and cosine (<£ - <p ) = cos n = - 1, 
as required for the maximum E A given by Equation 38. For any point in the F plane which does not 
fall on the line o x T , the cosine term will be greater than -1 and the value of e a will be less than 
that given by Equation 38. 

Examination of Equation 38 shows that the error E A increases as r and r x increase. There- 
fore, to determine the maximum value of E A as a function of r and r 9 only the maximum values of 
r and r x need be specified. Conversely, if e a is to be restricted to a value of less than or equal 
to a specified value, the maximum values of r and r must satisfy Equation 38 for that particular 
value of e a . 

In order to analyze the relation between maximum e a , r, and r x we must consider the inter- 
dependence of the maximum values of r and r x caused by the limitations of a finite lens aperture. 

In the discussion, it is assumed that diffraction effects at the end of the lens aperture are negligible. 
Figure 10 shows the extreme rays which can pass through a lens aperture of radius R L to reach 
the points at the distance r max from the optical axis. It should be apparent that any ray parallel to 
but above the upper extreme ray, or parallel to but below the lower extreme ray, will be outside 
the lens aperture and will not pass through the lens. Thus any signal point outside the ray defined 
by r lmax in Figure 10 cannot contribute to both of the points + r max and -r max . For example, if the 
signal area were extended upward beyond the r imax limit, the additional signal interval cannot con- 
tribute to the spectral point at +r max since the necessary light path will fall outside of the lens aper- 
ture. Under these circumstances the amplitude at the spectral point at +r max will not correspond 
to the entire signal but only to the interval below the +r lmax limit. From this example it is apparent 
that the r lmax limit given by Figure 10 defines the maximum signal interval over which every point 
contributes to the spectral points at ±r max ; 

The dashed lines in Figure 10 represent the extreme rays to a spectral point at a distance r 
that is less that r max . It is seen that the extreme rays for such a case define a maximum signal 
interval longer than that obtained for r . This means that the signal interval defined by r, 
increases as the spectral range of interest, defined by r max , decreases. Thus it is seen that for a 
given maximum frequency (i.e., r max ), the maximum value of the signal interval r lmax is limited by 
the lens aperture. 

To derive an expression defining the relation between r max and r imax , the geometry of Fig- 
ure 11 is used. This figure represents the upper extreme ray and the principal ray contributing 
to the spectral point at +r max . Since the extreme ray must be parallel to the principal ray, the 
angles 9 are equal and because of a similar-triangles relationship 

**max _ *"lmax _ 

f ~ a ~ a + d 
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Figure 11 “Geometry for relation between r and r, 

^ * ma x | max 


From these relations, two equations for a can be obtained: 


a 



a - 



- d . 


Since the right hand sides of these equations must be equal, then 


- t"i , max 

r imax ^ f 


A lens is usually specified by its F stop which is defined as 


F = 


2R, 


Dividing both sides of Equation 39 by f , 


^lmax _ d r max 

T~ ~ T ~ T T~ 


From Equation 40 it is found that R L /f = 1/2F; substituting into Equation 41 


(39) 


(40) 


(41) 


r, 

1 max 


1 d ^max 

2F “ T ~T~ 


(42) 
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It is obvious from Equation 42, as well as from Figure 10, that r imax cannot be greater than the 
lens aperture radius r l . Equation 42 defines the maximum allowed signal aperture radius r lmax 
because of the limitations of the lens aperture. In practice, the size of the signal aperture is 
specified by either physical consideration or a desired size format. The terms in Equation 42 can 
be rearranged to define the maximum spectral term r max as 


r 

max 


f 


1 r imax 

2F r~ 

d 

f 


(43) 


Equation 43 defines the maximum allowable r for a given r lmax as determined by the restriction of 
a lens aperture. By rearranging terms in Equation 38, a second expression can be obtained speci- 
fying the limitations on r max required for an allowed error E A : 


r 


max 


f 


1 r imax 

2 T~ 



(44) 


To demonstrate the application of Equations 43 and 44, consider the case for r lma /f = 1/5 
(e.g., for f = 100mm, r t = 20mm). Figure 12 is a plot of' r max /f as a function of d/f for the speci- 
fied input aperture, r x /f = 1/5. The curves labeled F = 1.4 and F = 2 correspond to Equa- 
tion 43 for the specified values of F. The 
curves labeled E A = 0.02, E A = 0.01, and e a 
= 0.005 correspond to Equation 44 for the 
specified values of e a . The f curves specify 
the upper limit on r/f due to the lens aperture, 
and the E A curves specify the upper limit for a 
given accuracy of the approximation. For a 
chosen value of d/f, the value of r max /f must be 
below both the F and E A curves which apply to 
the particular system being considered. 

For example, consider the case in which a 
lens with F = 2 is to be used and the maximum 
error to be allowed is E A = 0.02. The greatest 
value allowed for r/f corresponds to the point 
A which is the intersection of the F = 2 curve 
and the E A = 0.02 curve. The value of d/f 
= 0.5 would be selected to obtain the value 
r max/ f = 0.1. For any other value of d/f, the 
limit on r max /f would be less than the maximum Figure 12— Limitation on maximum spectral term r max /f. 
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at point A. For d/f less than 0.5, the E A = 0.02 curve specifies a tighter limit on r max /f, while for 
d/f greater than 0.5 the F = 2 curve limits r maj{ /f. Of course any combination of r max /f and d/f 
corresponding to a point below the curves is allowed; the curves only define the upper limit on 
r max/ f for a g iven value of d/f. 

The specification of a desired value of the maximum error e a is not readily determined in 
practice. Since the error in the contribution from (x^ y x ) varies from point to point, the total 
effect of the error cannot be determined unless the integration of the exact expression given by 
Equation 20 can be evaluated. To circumvent this difficulty, a more or less logical selection of 
parameters is considered, and the maximum E A specified by these parameters is determined. The 
value of e a determined will then specify the maximum error in the approximation be considered. 

In Figure 12, the dashed line corresponding to r max /f = 0.14 represents the limit determined 
in the last section for an accuracy of better than 1 percent in the linear relation between spectral 
frequency and back focal-plane coordinates. From the previous discussion of Equation 38 it is 
noted that the error E A decreases as d/f increases. Therefore, if it is not desired to lower the 
previous limit of r max /f = 0.14, the error, e a , can be improved (decreased) by using the largest 
possible value of d/f. Referring to Figure 12, note that for a lens with f = 2, the maximum value 
of d/f is 0.3, which is given by the intersection of the dashed line r/f = 0.14 and the f = 2 curve. 

As noted in Figure 12, the value of E A at this point is 0.037 (from Equation 38). If a lens with f = 1.4 
was considered, Figure 12 shows that d/f can be increased to a value of 1.1 and the error to e a 
= 0.023 can be reduced. Thus the usual result is obtained that the lens of lower F provides the 
better characteristics (lower F implies larger lens aperture for given focal length). In addition to 
having the larger error E A , the F = 2 lens restriction of d/f = 0.3 presents practical problems— the 
lens mount and input aperture mount must be designed to allow for a small spacing (d = 3 cm for 
f — 10 cm). 

In practice lower values of r lmax /f and r/f may be satisfactory. In such cases, the error e a 
would be less than the 0.023 determined here, and higher-F lenses (smaller lenses) may be used. 

An extreme case has been considered here and what amounts to an extreme error, e a = 0.023 
(or 2.3 percent), has been determined. This extreme value of the error, introduced by the Fourier 
transform approximation, is quite reasonable and should be sufficient justification for using the 
transform approximation in most applications. 

The following sections consider the phase of the light distribution in the back focal plane F' . 

It will be shown that selecting d/f = 1 has advantages in reducing the phase factors not associated 
with the Fourier transform. In Figure 12 the point corresponding to d/f = 1 and r max /f = 0.14 is 
shown to have an error value E A = 0.024. Thus it is seen that for the extreme case considered in 
the discussion above, reducing the value of d/f by one tenth increases the error by 0.001.- Such a 
slight increase in the error e a is quite reasonable in terms of the advantage gained in the phase 
approximation treated in the next section. In addition, since the location of a plane at a value of 
d/f can never be completely accurate, the displacement from the F = 1.4 curve allows a safety 
margin of +10 percent allowable error in the location specified by d/f - 1, without exceeding the 
limitations imposed by the lens aperture. 
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Thus it has been shown how Equations 38, 43, and 44 can be used to determine and/or specify 
the parameter limits and the accuracy of the Fourier transform representation: 


A(x, y) 


_ ie ikR(x,y) 

' X(f + d) 



y ) e -i 2 w (P x i + oyi) 


dXj dy j . 


(45) 


It has been shown in particular that for the maximum values r max /f = 0.14 and r lmax /f = 1/5, and 
the desirable choice of d/f = 1, the worse possible error in the amplitude values given by Equa- 
tion 38 is 2.4 percent. Since the term neglected in Equation 20 is negative, the approximate ampli- 
tude given by Equation 45 will be higher than the exact values by no more than 2.4 percent. 


FOURIER TRANSFORM REPRESENTATION OF OPTICAL IMAGING 

By limiting the area of consideration in the input plane, 


and in the back focal plane, 


/ r 1 max 1 \ 

\~f~ < 5/ ' 


( r max \ 

— < 0 . 14 ) 


it has been shown that the light amplitude distribution, A(x, y) , in the back focal plane of a lens is 
given with reasonable accuracy by 


-~i. e ikR ( x,y) 

A ( x > y) = \(f +d) F (P’ q) > 


(46) 


where 


and 


_ X _ y 

p - xf q ■ xf 


(47) 


f 2 + df + x 2 + y 2 

R(X ’ V) " ( f 2 + X 2 +y2)l/2 


f 2 + df + r 2 
( f 2 + r 2 ) 1/2 ’ 


(48) 


F(p. q) 



'i) 


^-iaw^pXj+qy,-) 


dXj dy j 


(49) 
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Figure 13— Two Sens optical imaging system. 


As given by Equation 49, f(p, q) is the two- 
dimensional Fourier transform of the light 
amplitude distribution in a plane 

perpendicular to the optical axis and at a dis- 
tance d in front of the lens. 

As pointed out in the discussion of Equa- 
tions 22, 23 and 24, the phase term e ikR ( x ’ y > is 
of concern only when a second lens is introduced 
to produce an image as shown in Figure 13.* 

To simplify our notation, k is introduced 
and defined as 


K 


. ie ikR(x,y) 

A.( f + d) 


(50) 


and Equation 46 is written 


A(x, y) - KF(p, q) . (51) 

In Figure 13, the light amplitude distribution A' (x t , yj in the input plane F is given as s(x l , y^. 
Using Equations 49, 50, and 51, the light amplitude distribution A(x, y) in the plane F' (back focal 
plane of lens L t ) is given as 


A(x, y) 



fl ) e ’ i27r ( px i +qy i) dXi dyi 


(52) 


where 


p 




(53) 


K 


1 


. ikR (x,y) 
“ie 1 

M f i +d i) ’ 


(54) 


and 


( x > y) 


f i 2 + d i f i + f2 

( fj 2+ r 2)i/2 


( 55 ) 


*Only conventional optical systems are considered here. If holographic techniques are considered, the phase factor in Equation 46 
would determine the form of the interference pattern produced by A(x, y) and a reference signal. 
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Similarly, A(x, y) is the input signal to the lens L 2 ; the light distribution, o(x 2 , y 2 ) , in the output 
plane F" (back focal plane l 2 ) can be written as 


°( x 2’ y 3 ) 



e ”i27T(P # x+q^y) dxdy ( 


where 


and 



ikR 2 ( x 2 , y 2 ) 

M f 2 + d 2 ) ’ 


r 2 (v y 3 ) 


f„ s + f,d„ + r„ 


(f 


2 W) 


2 \ 1/2 


(56) 


(57) 

(58) 


(59) 


Substituting Equation 52 into Equation 56, an expression can be obtained for the output image 
0(x 2 , y 2 ) in terms of the input image s(x 1 , y 1 ) : 


°( x 2- y 2 ) 


.a- 


i2 7T (P ; x +q f y ) 



”i2 7T(px j + qy i) 



(60) 


Assuming that the function s(x 1 , y x ) allows the order of integration to be reversed, then Equa- 
tion 60 can be rewritten as 


°( x 2 ’ y 3 ) 



vi) dx i dy i 



‘i2' 7T ^ P x !+qy j+p'x+q 


'y) 



( 61 ) 


The integral within the brackets is complicated by the presence of the factor K x , which contains an 
exponential dependent upon x and y . Limiting the values of x and y (i.e., r/f ) so that the phase 
variation in K t can be considered negligible, then the K x factor can be taken outside the integrals 
so that 


°( x 2- y 3 ) 



y i) dx i dy x 



(62) 
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Now consider the integral within the brackets and substitute for p, q, p' and q’ the values from 
Equations 53 and 57 so that 


j e -i2^px I+qyi+ p.x +q . y) dx dy 


I 


^i27T(X jA f i+x 2 Af 2 ) 


x dxJe' i27r ( yiAf > +y ^ A ^) y dx . 


( 63 ) 


Up to this point limits of integration have not been mentioned. Because of the presence of an 
aperture in an optical system, the signals exist only over a finite range of the aperture coordinates. 
However, since the signals are zero outside this range (e.g., s(x 1 , yj «= 0 outside the aperture in 
the F plane), the contribution to the integral beyond the aperture limits will also be zero. Thus 
the limits of integration are from to +co. These limits are in agreement with the Fourier trans- 
form integrals. 

Now, the Dirac delta function can be defined by the integral equation 

S(x-a) = | 


Comparing each of the integrals on the right side of Equation 63 with the integral in Equation 64, 
then 



-i2'n^px 1 +qy j + P' x + 


q,y ) dx dy 




= \ 2 f 2 


slx i + t: 






(65) 


The second step of Equation 65, uses the identity 


8 (ax) 



S(x) . 


Equation 65 is valid only if x and y range from -<» to +°°. hi optical systems this is not the case 
since the range of the coordinates x and y is limited as demonstrated in the previous discussion. 
However, Equation 65 is assumed valid to simplify the discussion. Substituting Equation 65 into 
Equation 62 gives 


°(v y a ) 


K 2 Kj a. 2 f 




'i) 


X ! + 


r a 


*1 


dx. 


dy x . 


( 66 ) 
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Now, the sifting property of the Dirac delta function which is defined by 


I 


F(x) 8(x + a) dx 


F(-a) , 


can be made use of. Applying this property of the delta function to Equation 66, it becomes 


0 (v y 2 ) 


k 2 k x X 2 S|- 



(67) 


In deriving Equation 67 it was assumed that K 1 was approximately constant The factor K 2 is 
variable only in phase, as seen by referring to Equation 58. Since only intensity is seen or meas- 
ured, the phase variations of k 2 can be ignored, and Equation 67 can be interpreted as giving the 
output image, o(x 2 , y 2 ) , in terms of a proportionality factor ( K 2 K 1 \ 2 f f ) which multiplies the orig- 
inal input signal, S, expressed in the new coordinates 

( f ! f l 

V ^ v ' r a y 2 



To clarify the significance of these new coordinates, consider the relation between the two sides 
of 


s ( x i’ 



( 68 ) 


Since the two sides of Equation 68 correspond point for point, then the ( x 1# y x ) and (x 2 , y 2 ) co- 
ordinates are related by 


i 1 1 

x i = ~t 2 x 2 and y i - " t; y 2 ■ (69) 

Equation 69 represents the fact that a signal point which was originally at the coordinates x l and 
y will be imaged to the point at 


X 


2 


77 x i and y 2 


f 2 

-77 y ! ■ 


The magnification in an optical image is defined as the ratio of the imaged coordinate of a point to 
the original coordinate such that 


m 


X 



f 


2 
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and 


^ 2 _ 2 

m y = 3^ = ~T;‘ (70) 

Equations 70 were written separately although it is apparent that here the magnification is the 
same in any direction. In some cases it is possible to obtain different magnifications in different 
directions (e.g., cylindrical lens system). Equations 67 and 70 show that the output image is pro- 
portional to the input image with a change in scale. The minus sign, which appears in Equation 67 
and 70, represents an inversion of the image. 

In many applications, there is no requirement for a magnified image. In such cases, lenses 
of equal focal length q = f 2 could be used, and a magnification of m = - 1 would be obtained. 

For the case f = f - f , Equation 67 becomes 

°(V Ya) = K 2 K t A . 2 f 2 s(-x 2 , -y 2 ). ( 71 ) 

Thus, for equal focal length lenses, the output image o(x 2 , y 2 )is proportional to an inverted replica 
of the input signal. 

It is for this case, f x = f 2 = f , that the optical imaging process can be described as consecu- 
tive Fourier transforms. This can be shown by replacing f by f in Equation 57 and substituting 
for p' and q' in Equation 56 to obtain 


0(x 


■) = 


= K 


JJ 


A(x, y) e 


- i 2 77 (X 2 xA f + y 2 y Af ) 


dx dy 


(72) 


By replacing f by f in Equation 53, then 


X 

p ~ JJ ' 


_ y 

q - Xf ’ 


dx - \f dp , 


dy ~ \f dq . 


(73) 


Substituting Equation 73 into Equation 72 


°( x 2’ y 2 ) 


K 2 X2 


f 2 JJ A(x, y) e" i2 ^ P ^ +qy ^dpdq 


(74) 
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Using Equations 51, 52, 71 and 74, the two step process of optical imaging can be expressed as 


A(x, y) - KjFCp, q) 



y i) 


-i 2 ^(px l+qyi ) ^ 


l dy X 


( 75 ) 


and 


°( x 2> y a ) = K 2 Kj X 2 f 2 S(~x 2 , -y 2 ) 


k 2 K^ 2 f 2 


II 


F(p> q) e 


"i2 7T (px 2 + qy 2 ) 


dp dq 


(76) 


where 


X 

kf ' 


_ y 
~ XT 


. ikRj (x,y) 

~ie 1 

( i ~ ^(f + d d 


_ ie ikR 2 ( x 2-y2) 


= df+d 2 ) ’ 


f 2 + fd. + 


f 2 + fd 2 + r 2 2 


R 1 ( X ’ y ) (f2 +I . 2 )l/2 ’ R 2 ( x 2 ’ y 2 ) ^ f 2 + T 2 ^ 


1/2 


r 2 - x 2 + y 2 > 


and 


^ 2 f y 2 
2 y 2 


The last expression in Equation 76 assumes that the factor k x can be considered constant in phase 
over the range of the values p and q (i.e., x and y). This approximation is the subject about to be 
considered; shown here are the advantages of the resulting expression. To appreciate the sig- 
nificance of Equations 75 and 76, consider the standard Fourier transform equations in the nota- 
tion used here : 


F(p. q) - j 

js(x 1 , yi) e- i 2 ^^ l )dx 1 dy 1 , 

(77) 

S K> y i) = 

JJ F(P , q) e +i277 ( PX ‘ +qyi) dpdq , 

(78) 

S(-x x) -y x ) 

= J"J"f( P , q) e ~ i27T ( px i +qy i) dp dq . 

(79) 
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Equation 77 is usually referred to as the Fourier transform while Equation 78 is the inverse Fourier 
transform. Note that the exponent of Equation 78 is positive and that of Equation 79 is negative. 

Since the optical transform produced by a lens has a negative exponent, the inverse transform de- 
fined by Equation 78 never appears in optical systems. The second lens in an optical system, such 
as that of Figure 13, produces a Fourier transform of a Fourier transform as represented by Equa- 
tion 79. Note that the inversion, or change of sign of the coordinate, is introduced by the second 
Fourier transform, whereas an inverse transform would not invert the signal. Thus, comparing 
Equations 75 and 76 with Equations 77 and 79, it is noted that the optical imaging process of two 
lenses is described by two successive Fourier transforms relations. Except for determining the 
absolute amplitudes involved, the constants in front of the integrals of Equations 75 and 76 do not 
affect the form of the variations. In most cases only the relative amplitudes are of interest, and 
the constants are dropped. 

The advantage of the Fourier transform representation described by Equations 75 and 76 can 
be shown by considering the introduction of a filter in the F' plane. If the transmission character- 
istics of the filter are known, then a function M(p, q) can be determined which represents the frac- 
tion of incident light amplitude passed at each coordinate corresponding to the values of p and q. 

The filtered output 0 f (x 2 , y 2 ) is then given by Equation 76, if F(p, q) is replaced by M(p, q) F(p, q) 

Of ( ^ 2 * y 2 ) = K 2 K l X2 f 2 JJ" M ( p> q ) F ( p ’ q ) e l27r ( pX2 + cly 2) dpdq • (80) 

Thus the specification of a filter for a particular application can be determined uniquely when the 
Fourier transform representation is used. 

It is noted at this point that the Fourier transform representation of Equations 75 and 76 re- 
quire the use of lenses of equal focal length. The specification of a filter for the case of unequal 
focal lengths is exactly the same; however, the Fourier transform relation of Equation 76 is modi- 
fied by introducing a factor of f 1/ /f 2 in the exponent to account for the magnification. 

This additional f 2 /f 2 results in a magnified filtered image which contains the same information 
as the filtered image, o f (x 2 , y 2 ) , given by Equation 80. The only difference (neglecting the phase 
of k 2 ) is in the scale. Throughout the remaining part of this report, the special case of equal focal 
length lenses is considered to simplify the analysis. 


ELIMINATION OF UNDESIRABLE PHASE VARIATIONS 

Having now seen the significance of the Fourier transform representation of optical imaging, 
it is possible to consider the approximation involved in the derivation of Equation 76. The actual 
relation corresponding to Equation 76 can be written as 


0(x 


.) = V 2 


J K ‘ 


F(P> q) e 


-i277(px 


f qy dp dq . 


(81) 
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The term K t appearing in the integral was defined as 


where 


ikRjCx.y) 

K i = \(f + dj ’ 


( 82 ) 


R i ( x - y) 


f 2 + fdj + r 2 
(f 2 tr 2 ) 1 / 2 


(r 2 = x 2 + y 2 ) • 


(83) 


Comparing Equations 81 and 76, it is apparent that Equation 76 is valid only if the phase variation 
of K : can be neglected over the range of the integration variables p and q (or x and y). This re- 
quirement is now analysed. 

For reasons that will become evident, R (x, y) is expressed in the form 

R t (x, y) = (f +d x ) + P(f-d t ) + Qf , 

= f(l +P+Q) + d t (1-P) . (84) 


Equation 83 can be rewritten as 


_ f 2 + r 2 d i f 

R ! ( x - y) - ~ - + 


(f 2 + r 2 ) 1/2 (f 2 + r 2 ) 1/2 




• 1/2 


(85) 


Equating the coefficients of f and d in the final forms of Equations 84 and 85 ? then 


and 


Solving for p and q 


1 - P 


/ r 2\" 1/2 
\ 1+ f 2 ) 


1 + P + Q 


/ r 2\ 1/2 

\ 1+ f 2 j 


( * 2 \ 1/2 

l 1+ f 2 J 
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Q 


2 


/ r 2 N 1/2 / r 2\~ 1/2 

V + f 2 ) + l 1 + f 2 ) ' 2 


1 + 


2f 2 


( r 2V 

y f 2 ) 


1/2 


Substituting for p and q, Equation 84 can be written 


(x, y) = (f + d 1 ) + (f-dj 


(*^r 


+ 2f 


1 + 


2f 2 


/ 2 \ 1/2 

K) 


- 1 


(86) 


The K x term given by Equation 82 can be rewritten using Equation 86 so that 


K 


i 


: ie ik ( f+d i>' 

.H {+ y. 


ik(f-d 1 )[l-( l+r 2 /f 2 )~ 1/2 ] 


i2kf [( l+r2/2f 2 )( l+r 2 /f 2 )~ 1/2 -l] . 

e 


(87) 


The terms grouped within the first brackets of Equation 87 are constant and therefore can be taken 
out from under the integral sign in Equation 81. The remaining exponentials in Equation 87 are 
phase factors which depend on the variables of integration. The exponential of the remaining terms 
must be limited so that the phase variations can be considered negligible. 

Since restrictions must be considered on the value of r/f , so that the phase terms can be con- 
sidered negligible, then the exponentials of Equation 87 must be simplified by expanding in power 
series of (r 2 /f 2 ) while the first terms of the expansions are all dropped. Expanding the exponent 
of the first exponential gives 


ik(f-d 1 )[l-(l + r 2/f 2 ) _1/2 ] _ e ik(f-d 1 )[l-{l-(r 2 /2f 2 )+3/8(r2/f2) 2 - • •}] . 

-v ik( f-d , ) ( r/f ) 2 / 2 

~ e 

Similarly, by expanding the second exponential 

e i2kf[(l+r 2 /2f 2 )( l+r2/f2)-l'' 2 -i] _ i2k{ [~l+( 1 + r 2 /2f 2 ){l-( r 2/2f 2 )+3/8( r 2/f 2 ) 2 - • •}] . 

~ e ikf(r/f) 4 /4 . 

Substituting these approximate terms in Equation 87 gives 

_ -j e lk ( f d l) ik(f“dj)(r/f ) 2 / 2 ikf(r/f) 4 A 
K i - ^(f+dj) e e 


( 88 ) 
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For values of (r/f) less than our previous limit of 0.14, the. approximations in each of the phase 
terms is accurate to within 2 percent of its exact values. It should be noted that neglecting terms 
in the exponentials is valid only since phase variations less than one cycle are to be considered. 

Since the desire is to eliminate the phase variations of K v Equation 88 indicates that the opti- 
mum choice of the distance d 2 equal to f eliminates the first phase term. Thus for the case when 
d i is chosen equal to f , Equation 88 can be reduced to 


K i = 


■ie i2kf 

2^f 


»ikf(r/f ) 4 / 4 


(for dj = f ) 


(89) 


If Equation 89 is substituted for K in Equation 81, then 

-iK 2 A.fe i2kf 

°( x 2>y 2 ) = 2 — 


e 


ikf(r/f) 4 /4 


F(p- q)e 


i 2^(p x 2 +qy 2 ) 


dpdq 


(90) 


When the value of (r/f) is limited so that the phase factor appearing under the integral of Equa- 
tion 90 can be neglected, the output image o(x 2 , y 2 )is given by 


°( x 2-y 2 ) 


K x K 2 X 2 f 2 


II 


F(p, q) e 


-i2w(px 2 +qy 2 ) 


4) dq 


(91) 


where 


”ie l2kf 

K i = 2 x 1 ■ 


Equation 91 is identical to Equation 76 which has been shown to be the desired form for the Fourier 
transform representation of optical imaging. 

To derive a specification for the maximum limit on (r/f) which allows the variable term in k 
to be neglected, consider the effect of the phase term for a particular F(p, q) : 


F(p, q) 


S(q) 


A 0 s (p) + ?{ 8 (p-p 0 ) + 8 ( p + p o)} 


(92) 


The locations of the frequency terms contained in Equation 92 are diagrammed in Figure 14. Since, 
by definition, the delta function S(q) is equal to zero for q unequal to zero, Equation 92 represents 
the spectrum of a signal which varies in only one dimension. That is, there are no frequency com- 
ponents in the y direction; therefore, the signal is constant with respect to the y coordinate. Rather 
than interpret Equation 92 as the spectrum of a particular signal, it can also be assumed that only 
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Figure 14— Location of frequency terms 
In the spectrum plane. 


three sample points of a more general spec- 
trum are being considered. Since there is 
nothing to single out the x direction in an op- 
tical system, the analysis will apply to a set 
of spectral points along any radial axis in the 
frequency plane as indicated by the r axis in 
Figure 14. This is obvious if it is considered 
that arbitrarily selection of any orientation can 
be made for the x- and y- coordinate axes. It 
can also be shown that the terms used to specify 
a maximum limit on (r/f) also apply to the 
general case. Thus the results are simply 
interpreted for the special case of Equation 92 
as a general criteria for neglecting the unde- 
sired phase factor in Equation 90. 


Substituting for F(p, q) as given by Equation 92 and applying the interpretation of the discus- 
sion above, Equations 90 and 91 can be simplified to 


°K) = Kje^^^0 4 A[ AoS(p) + || s ( p _ p J + s ( p + p jj 


~i2^px 

e * dp 


( x 2 ) - kJ A 0 S(p)+ \ {s (p -p 0 ) +S(p +p 0 )| e l27rpX2 dp 


(93) 


(94) 


The new factor K in Equations 93 and 94 is defined as 

-iK„ \ f e i2kf 

TF ~ 


(95) 


The integral with respect to q was taken by applying the sifting property of the delta function S(q). 
The y 2 dependence has been dropped on the left side of Equations 93 and 94 since the output image 
varies only with respect to the x 2 coordinate. The left side of Equation 93 is underlined so that 
we can identify the ensuing results of Equations 93 and 94 throughout the remainder of our 
discussion. 

The first exponential in Equation 93 can be rewritten using the definition p = x/\f , while the 
integration with respect to p is performed simply by applying the sifting property of the delta 
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function: 


>(x 2 ) = K je ikf ^ p ) 4 A[A 0 5(p) + | {s(p+p 0 )+8(p-p 0 )} 


-i27Tpx r 


dp 


- t/Ta . B ikf(X.p 0 ) 4 A f i2wp 0 x 

- K |A 0 + ~2 e \e 2 + e 


i27Tp o x 2}j 


°(x 2 ) " k[a q +B e lkf<A ' P0) ^cos 2 7rp 0 x 2 ] • 


(96) 


The result after integrating Equation 94 is similar to Equation 96 except for the exponential which 
appears in Equation 96: 


°( x 2 ) = K[A 0 +Bcos2tt Po x 2 ] . (97) 

Equations 96 and 97 represent the amplitude distribution of the output image produced by the 
frequency-plane distribution F(p, q) given by Equation 92. Equation 96 represent the output o(x 2 ) 
when the phase factor is considered and Equation 97 represent the output o(x 2 ) when the phase 
factor is neglected. Comparing Equation 96 and 97, a criteria for neglecting the phase factor is 
still not very apparent since the significance of the exponential is not very clear. 

If the observation of the output image is considered, then the intensity, rather than the ampli- 
tude as given by Equations 96 and 97, must be dealt with. The intensities are given by the relations 


X ( x 2 ) = °( x 2 ) °*( x 2 ) = l °( x 2 )| 2 ’ 


(98) 


and 


l(x 2 ) = 0(x 2 )o*(x 2 ) = |o(x 2 )| 2 . ( 99 ) 

The starred terms in Equation 98 and 99 represent the complex conjugate of the unstarred terms. 
Using Equations 96 and 97 in Equations 98 and 99 respectively 

r kf "i 

I(x 2 ) = |k| 2 [A 0 2 +B 2 cos 2 p 0 x 2 + 2A 0 Bcos (X.p 0 ) 4 cos 2 tt p 0 x 2 J , (100) 
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and 


I(x 2 ) - |k|2[a o 2 +B2 cos2277 Po x 2 +2A 0 Bcos2^p 0 x 2 ] . ( 101 ) 

Now comparing Equations 100 and 101, it is found that the phase factor introduces a cosine factor 
which attenuates the cos 2 77 Pq x 2 component in the observed image. In the general case, the ampli- 
tude B of any one component will be considerably less than the component A 0 . Therefore, the third 
term in Equations 100 and 101 represent the larger of the two x 2 -dependent terms in the image in- 
tensity. It is desirable to limit the maximum value of p Q to obtain a value of cos (kf/4) (\p 0 ) 4 as 
near to one as possible so that the Fourier-transform representation used in deriving Equation 101 
can be considered a good approximation. 

The cosine term can be expressed as a function of x by the definition p = x/\f, and since the 
results will apply to the general case, x can be replaced by the more general notation r. Thus the 
cosine term can be expressed 


kf 


i^PoY 


kf 


(f ) 4 


( 102 ) 


In Equation 102 the frequency p 0 is given the general interpretation of a spatial frequency in the 
direction of an r axis (see Figure 14), and p 0 and r are related by p Q = r/\f . That is, Equation 102 
applies to the general case of a spectrum along any radial axis r in the back focal plane F\ The 
term (r/f ) can be expressed as a multiple of (4/kf )^ 4 by defining a factor m by the relation 


r 

T 



(103) 


Substituting for (r/f) in Equation 102, 


COS 




cos m 4 . 


(104) 


The ultimate limit on m can be recognized by noting that for m = ( 77 / 2 ) 1/4 the cosine term as given 
by Equation 104 is zero (i.e., cos 77/2 = 0). For this value of m the Fourier transform result of 
Equation 101 is completely in error with respect to the third term, since the cosine term present 
in Equation 100 is zero and the third term is eliminated. Thus for m = ( 77 / 2 ) 1/4 , the ideal two step 
Fourier transform representation yields a term which does not exist in the actual image given by 
Equation 100.* For values of m less than ( 77 / 2 ) 1/4 , the cosine of Equation 104 has nonzero values 


♦This accounts for the frequent neglecting of phase terms, only when less than 77/2 which appears in many references dealing with ap- 
proximate solutions of diffraction problems. 


36 



as shown by the curve of Figure 15. Selecting 
a limit on m is based on specifying how accu- 
rately the third term of Equations 100 and 101 
should agree. A value of m = 0 is necessary 
to have complete agreement between Equa- 
tions 100 and 101; however, m = 0 corresponds 
to a frequency p 0 = 0 which corresponds to a 
dc term only. Thus a compromise limit must 
be established between the limits m = 0 and 
m = (77/2) 1/4 . 

To determine the limitation onm,it is nec- 
essary to specify the desired accuracy of the 
third term in Equation 101 as compared to the 
third term in Equation 100. Again, the accu- 
racy of the approximation can be given in term; 

1 - cos m 4 1 

^ cosm 4 cosm 4 * (105) 

The percent error in the third term of Equation 101 is then +100 E^ percent compared to the exact 
term in Equation 100. Note that the error does not represent a fraction of the total image in- 
tensity. The error E^ corresponds only to a particular term in the image intensity. In the general 
case there would be a series of such terms and the maximum E 0 would be determined by Equa- 
tion 106. This E 0 would represent the maximum error in terms of the form of the third term in 
Equation 101 and would correspond to the term involving the highest frequency of interest. 

Specifying the maximum allowable error is rather arbitrary and will usually depend on the 
particular application considered. However, for an example a limit of 2% accuracy can be specified 
for the approximation, i.e., E^ = 0.02. Substituting in Equation 105 the condition for the maximum 
value of m is obtained: 



Figure 15— Graph of cos m 4 vs. m. 


of a fractional error E^ defined as 


n 4 * 


1.02 


= 0.98 


(106) 


From Figure 15 it is found that Equation 106 requires values of m less than or equal to 0.67. Using 
the maximum value m = 0.67, Equation 103 becomes 


(t) max = 0-67 (u) ■ (107) 

To obtain a numerical result for comparison with the previous limit, r/f ma)i = 0.14, a wavelength 
\ = 5461 x 10 “ 8 cm and focal length f = 10 cm is considered. Substituting in Equation 107 it is 
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found that 


(t) 

\t/ max 


( 4X\ l/4 
0-67 fen) 


0.67 


4x 5461 x 10 ~ 8 V /4 

2 rrx 10 / 


(tLx - 0-03 . 


(108) 


Equation 108 specifies the aperture limit in the frequency plane F' to assure an error limit of less 
than 2 percent caused by neglecting phase variations. Notice that this phase limit restricts the fre- 
quency aperture to approximately one fifth of the previous value of r/f max = 0.14, which was sufficient 
for the linearization and amplitude approximations. The corresponding frequency limit is given as 


1 ( r \ ^ 0.03 

Pmax = X \*F /max " 5461 x 1Q - 8 = 550 Cydes/cm . ( 109 ) 

Thus for spatial frequencies less than 550 cycles/ cm, neglecting the phase term to obtain the 
Fourier transform representation of Equation 91 introduces an error of no more than 2 percent in 
terms of the form of the third term in the image intensity of Equation 101. It is again pointed out 
that for most practical cases, the frequency capability of present input techniques restricts the 
possible frequencies to a lower value than that specified by Equation 109. 

It has been shown how Equations 103 and 105 are used to determine the error E^ for any fre- 
quency plane aperture with a radius defined by ( r max /f). It has further been shown for a particular 
case (k = 5461 A and f = 10 cm) that the limit, (r/f ) max = 0.03, provides an accuracy within 2 per- 
cent for the terms in which the phase variation appears. It has also been noted that this phase ap- 
proximation requires a tighter restriction on the maximum frequency terms. In fact, for examples 
used, the maximum frequency is one fifth of that allowed for an accurate amplitude approximation. 
Of course, this further restriction of the frequency range of interest will also improve the accuracy 
of the amplitude approximation. 

Figure 12 can be used to consider the amplitude error e a for the values d/f = 1 and r max /f 
= 0.03 assuming r x /f = 1/5. The point corresponding to d/f = 1 and r max /f = 0.03 is located 
below the curve corresponding to E A = 0.005. Therefore, the further restriction on r max /f required 
for the phase approximation reduces the error in the amplitude approximation to a value less than 
0.5 percent. This result shows that the restriction considered in this section not only provides a 
Fourier transform relation, which is accurate in phase, but also improves the accuracy of the am- 
plitude approximations previously considered. 

Within the limits presented in this section, the two- lens optical imaging process can be de- 
scribed by the equations 


A(x, y) 


hie 


i2kf 
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F(P> q) 
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-i27r(px 1 +qy 


dx 


i d yi 


(no) 
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and 


/ . ikR, [ ~ - f e l2kf ~ | ff -i2w(px ,+qy,) , , 

°(x 2 ,y 2 ) = e L 2 ( f +d 2 )J JJ (p ’ q)e dp dq (Hi) 

In practice only the variations in amplitude are of interest and the constant factors within the brac- 
kets are dropped, so that 


A(x, y) 


F(p, q) 



yj 


^“i27T(px 1 +qy 1 ) 


dx 1 dy t , 


( 112 ) 


°(x 2 , y 2 ) 


s(-x 2 , -y 2 ) 


ikR , 


Jf 


F(p, q) e 


-i27T(px 2 + qy 2 ) 


dp dq - 


(113) 


Equations 112 and 113 represent the form of the optical Fourier transform representation com- 
monly used. These equations describe the relative amplitude and phase variations of spectrum 
A(x, y) and image 0(x 2 , y 2 ) . Note that the phase term e lkR 2 is retained in Equation 113. This 
factor has no effect on the image intensity since multiplication by the complex conjugate eliminates 
this term. However, if the image o(x 2 , y 2 )is to be processed further by another lens, the effect 
of the phase factor e lkR 2 must be considered. In such cases, the criterion for neglecting the vari- 

ikR 

ation in phase due to the factor e 1 must also be reevaluated since the criterion developed above 
was based on image intensity effects. 


OPTICAL CORRELATOR SYSTEMS 

Considered now is a three lens optical system as shown in Figure 16. In this system the sig- 
nal plane F is assumed to be in the front focal plane of lens l , and each of the other lenses (l 2 
and l 3 ) is located so that its front focal plane coincides with the back focal plane of the preceding 
lens. With this configuration the amplitude distribution corresponding to the input signal to each 
lens is the output signal in the back focal plane 
of the preceding lens and is a focal length in 
front of the lens. This location of the signal 
planes provides the advantage of eliminating 
the phase terms dependent upon the distance 
from the lens to the input plane as discussed in 
relation to Equations 88 and 89. The optical 
system in Figure 16 consists of a two-lens im- 
aging system as discussed in the preceding 
section followed by a third lens that produces 
a Fourier transform of the light amplitude 


s(*i/ yd 



Figure 16— Optical correlator system. 
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distribution of the image. As pointed out at the end of the last section, the processing of the image 
°(x 2 . y 2 )by an additional lens involves its amplitude rather than the intensity; therefore, the phase 
effects of each lens will be considered. Throughout this section it is assumed that the aperture 
limitations are sufficiently restrictive so that the linear frequency and amplitude approximations 
developed earlier are valid. The focal lengths of the three lenses are assumed equal to sim- 
plify the analysis; in the general case, unequal focal lengths would introduce magnification or 
demagnification. 

In the optical system of Figure 16, lens L t produces a light amplitude distribution in its back 
focal plane F' which is proportional to the Fourier transform of the input signal, s(x 1 , y l ) , except 
for a multiplicative phase factor. As discussed in reference to Equations 112 and 113, all constant 
factors are dropped and only terms which vary with respect to the coordinates in the four signal 
planes of interest (f, f' , F", f 0 ) are retained. Using only the variable exponential in Equation 89 
for Kj , the amplitude in the F' plane is given by Equation 75 which can be written 

A(x, y) = e lkt(r/,) 4 / 4 F(p, q) = e lkf( r/f > / 4 JJs(x r yj e 12 1 9y ^ dx x dy x . (114) 


For our development of a correlator it is advantageous to introduce notation for the signal s(x l , yj 
which accounts for displacement of the signal from some reference position. Referring to Figure 
17, signal point A can be considered at the new position A' . If A is a point of the signal s(x 1 , y t ) , 
the light amplitude at A is s(x A , y A ). Since A' is the same signal point as A (it has only been moved), 
the light amplitude at A' must also be s(x A , y A ). The coordinates of the point A' are x x = x A +r x 
and y j = y A +r y . Thus the notation for the signal must be such that if the coordinates x x , y x are 
substituted for the point A', then s(x A , y A ) is obtained. The required notation is s(x 1 ~ T X > y 1 ~ T y ) 
as can be seen by substituting the values of x 1 and y 1 for each of the points A and A' . In either 
case the signal amplitude is s(x A , y A ). Using this new notation for a signal, Equation 114 can be 
rewritten as 



Figure 17— Displacement of a signal point 
in the input plane F. 


A(x, y) = e ikf ( t/f ) 4 Afi(p | q ) 

(115) 


= e ikf(r/fr 


The displacements r x and r y are positive when 
the displacement is in both the positive x x and 
positive y x directions. 

In Equation 115 the function F(p, q) repre- 
sents the Fourier transform of the displaced 
function s(x i -r x> y x -r ) . Thus the f(p, q) in 
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Equation 115 includes all the information regarding the signal including its displacement. From 
Fourier transform theory, the transform corresponding to a displaced signal such as in Equa- 
tion 115 differs from the transform F(p, q) of the undisplaced signal of Equation 114 by an exponen- 
tial phase term, e 127T < V * P+T y q ) . This principle need not be of any further concern; it is pointed out 
only to emphasize that the F(p, q) in Equation 115 corresponds to the displaced signal s(x ; -r x , 

yi" T y)- 

The amplitude distribution given by Equation 115 appears in the plane f' and represents the 
input signal to the lens l 2 . Lens L 2 performs a Fourier transform operation on A(x, y), and the 
image amplitude in the plane F" is given by 

0( x 2’ y 2 ) = ^( X 2 ’ y 2 )jj elkf<r/f)4//4 F (P’ 3) e ^ 2 9y2 ’ > d P d< l • (116) 


Equation 116 corresponds to Equation 81 except that only the variable terms of Kj and k 2 have 
been retained. The exponential appearing in the integrand corresponds to the variable term in Kj 
as discussed above. The function <f>(x 2 , y 2 )in front of the integral represents the variable part of 
K. . From the definition of K„ given under Equation 76 the variable part of K. is obtained from the 

2 ikR * ^ 

term e 1 2 where 


f2 + fd 2 + r 2 2 

■p z: - — — , and r 2 — 

K 2 ^f2 + r 2^1/2 n 2 X 2 Y 2 


Since R 2 has the same form as r i5 r 2 can be expanded in the form of Equation 86 


r 2 ( x 2>y 2 ) - ( f + d 2 ) + ( f _d 2 ) 


i- 


i + 


r 2 2\V2 


+ 2f 


1 + 


2f 2 


1 + 


r 2 ' 2 

T 2 


- 1 


(117) 


(118) 


The first term and the -1 term in the last bracket of Equation 118 are constant and can be dropped 
since only the variable part is of interest. The second term vanishes since d 2 = f in the system 
being considered. Thus the only variable term in Equation 118 is the fraction in the brackets of 
the third term. The function d>(x 2 , y 2 ) is therefore given by 

^ = e i2kf t(l + r 2 2 /2f 2 )(l+r 2 2 /f 2 r 1/2 ] . 

The variable part of K 2 given by Equation 119 was derived from the complete expansion of 
rather than from an approximate expansion analogous to Equations 88 and 89, since the aperture 
restrictions necessary for the validity of Equations 88 or 89 would require a signal and image 
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aperture much smaller than that normally desired in optical systems. For an image aperture de- 
fined by (r 2 /f) max = 0.14, and f = 10 cm, and X. = 5461 X 10 ~ 8 cm, the phase term y 2 ) can 
introduce phase shifts as great as 38 n radians (19 cycles). It was pointed out that the phase ap- 
proximation of Equation 89 was accurate with 2 percent. For the image aperture considered here, 
this phase inaccuracy can be of the order of 0.4 cycles. This magnitude of phase error may not be 
negligible, and, therefore, the more complete exponential was used in defining <£(x 2 , y 2 ) by Equa- 
tion 119. 

Returning to the image amplitude distribution o(x , y 2 ) given by Equation 116, the notation 
will be changed to take into account the possibility of image displacement corresponding to the 
signal displacement previously considered. In the last two sections it was pointed out that the 
imaged amplitude, o(x 2 , y 2 ), corresponds to an inverted replica of the input signal s(x x , y x ). This 
inverted property of the image applies to the image motion as well. That is, if the signal is dis- 
placed in the positive x x and y x direction, the image is displaced in the negative x 2 and y 2 direc- 
tions. Thus, if o(x 2 , y 2 ) corresponds to the inverted image of s(x x , yj, the displaced image cor- 
responding to s(x x - r x> y x - r y ) is obtained simply by reversing the sign of the displacement to 
obtain o(x 2 + r x , y 2 + r y ). Using the displacement notation for the image amplitude distribution, 
Equation 116 can be rewritten as 
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( 120 ) 


The final lens, l 3 , in Figure 16 operates on the light amplitude distribution appearing in its 
input plane f". For an optical correlator operation, a reference signal r(x 2 , y 2 ) is inserted into 
the plane f" in the form of an amplitude transmission function of a photographic transparency. In 
this case the light amplitude distribution operated on by lens l 3 is that which appears on the output 
side of the reference transparency. This light amplitude is given by the product of the incident 
light amplitude, o(x 2 + r x , y 2 '+ r y ) , and the reference transmission function r(x 2 , y 2 ). Thus the 
light amplitude distribution w in the output plane F 0 is given by 
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Since the system being considered terminates at the F 0 plane, the intensity will be detected, meas- 
ured, or recorded in the F Q plane. The intensity in the output plane is given by the product of 
Equation 121 and its complex conjugate. The complex conjugate product of the exponential in front 
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of the integral results in the cancellation of the exponential. Thus, the exponential in Equation 121 
can be dropped since it will not affect the detected intensity output. Equation 121, therefore, can 
be simplified to 
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Finally, if only the point located at the intersection of the optical axis with the plane f q (back focal 
point of l 3 ) is considered, s = t = 0 (i.e., x 3 = y 3 = 0), and Equation 122 reduces to 

W ( T x> T y) = JJ R ( X 2 > y 2 )°( X 2 +T X’ y 2 + T y) dx 2 d y 2 > (123) 

where 


W(t x , T y ) = w(x 3 = 0, y 3 = 0, T x , r y ) . 

Equation 123 corresponds to a two-dimensional correlation function that implies that the light 
amplitude at the back focal point ( x 3 = y 3 = 0) of the lens l 3 is given by the cross correlation of 
the reference, R(x 2 , y 2 ) , and the image amplitude, o(x 2 , y 2 ). Thus as the input signal is displaced, 
the variation of the light amplitude w (r x , r ) corresponds to the variation of the correlation func- 
tion with respect to the displacements r x and r . Note that the correlation function defined by 
Equation 123 involves the image amplitude o(x 2 , y 2 ) which is inverted with respect to the input sig- 
nal s(x x , y t ) . Therefore, if R(x 2 , y 2 ) is not a symmetrical function, it must be oriented correctly 
with respect to the image o(x 2 , y 2 ) rather than with respect to the input signal s(x , y ) . 

Considered briefly now is the implication of the steps from Equation 122 to Equation 123. This 
step in the derivation was accomplished by stating that only the single point in the output plane f o 
which lies on the optical axis (i.e., x 3 = y 3 = 0) could be considered. In practice it is physically 
impossible to isolate a single point. The best attempt that can be made is to restrict the light 
measurement or detection to a small area about the selected point. The light amplitude at points 
within this area (except for the one point on the optical axis) is given by Equation 122 rather than 
by Equation 123. The light amplitude distribution will not be uniform over the finite area of meas- 
urement because of the phase variation involved in the integral of Equation 122. For example, 
using a pinhole aperture 10 microns in diameter to define the detection area, the phase term in 
Equation 122 can vary as much as 4w radians (2 cycles) over the range of the image aperture (as- 
suming r 2max ~ 0.14f, f = 10 cm, K = 5461 x 10‘ 8 cm). The effects of the phase term in Equa- 
tion 122 is to reduce the light amplitude at points off axis since the contributions to the integral are 
not in phase. Therefore, the actual light available through a pinhole aperture located in the F 0 plane 
at x 3 = y 3 = 0 will be less than that found by assuming the light amplitude given by Equation 123 
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appears at all points within the pinhole aperture. This problem will not be considered any further 
here since the analysis would depend on the type of photodetector or measurement technique used. 

It will be assumed that the variations involved are small enough so that any measurement will yield 
values proportional to the square of the amplitude given by Equation 123. 

As pointed out above, the correlation function defined by Equation 123 involves the image am- 
plitude o(x 2 , y 2 ) rather than the single amplitude S(x 1 , yj. As defined by Equation 120 the image 
amplitude contains phase terms not present in the signal. A correlation operation can be per- 
formed based on the image as given by Equation 123; however, the reference signal r(x 2 , y 2 ) 
would have to be selected in terms of the image o(x 2 , y 2 ) including the phase terms. The correla- 
tion function obtained would correspond to a distorted signal rather than the actual signal s(x 1 , y t ) . 
The presence of distortion caused by the phase terms in Equation 120, therefore, complicates the 
analysis and determination of the correlation process. For example, the image amplitude o(x 2 , y 2 ) 
will be complex (phase variation as well as amplitude), and for complete correlation a complex 
reference signal is required. Such reference transparencies are difficult to produce. The phase 
distortions are commonly neglected and a reference signal is selected on the basis of an ideal 
image (no distortion) of the input signal. Analysis is now made of such a system to determine the 
effects of the undesirable phase terms present in our equations. 

Consider a signal which would produce an ideal image amplitude defined by 
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Equation 124 defines a signal image composed of a series of cosine harmonics in one dimension. 

A one -dimensional signal has been chosen to simplify the analysis. Referring to Equation 120 it is 
found that each frequency term in the image has a phase term e lkf ( Pn ) 4 associated with it, and 
the image also has a phase term <£(x 2 , y 2 ) associated with it. Thus, for the actual image, 

°( X 2 +T x) = ^( X 2 )^] B n e lkf<XPn) /4c Os 277 p n (x 2 +T x ) . 


A reference signal without phase can be considered such that 
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The reference signal r(x 2 ) defined by Equation 126 has been selected to have the same cosine 
harmonics (m = n) as the imaged signal being considered. Note that the reference signal defined 
by Equation 126 does not contain the phase terms present in Equation 125. The product of reference 
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and image for the ideal image of Equation 124 is given by 
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The product of reference and image for the actual image of Equation 125 is given by 
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Substituting Equations 127 and 128 into Equation 123, the result for the correlation function of the 
ideal image is 
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and for the correlation functions of the actual image: 
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Comparing Equations 129 and 130 it is found that the term e lkf ^ p "^ ^ affects the phase of each 
term in the double summation. From the previous discussion of frequency limitations with respect 
to this phase term it can be seen that applying the limitation developed for imaged intensity limits 
the phase variation of this term to approximately 7 degrees, h l summing terms which are not in 
phase, the result will be less than summing the same terms in amplitude only. Thus the presence 
of the phase term e lkf ^ p "^ /4 has the effect of reducing the value of w (r x ) in Equation 130 as com- 
pared to Equation 129. However, since the maximum phase will be about 7 degrees, the difference 
due to this term will be small. The phase term <£(x 2 ) appears in the integral of each term in the 
sum and has the same effect on the integral (can be considered as summation) as the phase term 
discussed above had on the summation. However, as discussed above, the phase variations of <i(x 2 ) 
ranges over 19 cycles and the effect on the value of tbe integral will be correspondingly greater. 

The actual magnitude of the reduction in w(r x ) because of these phase terms is difficult to evaluate 
in general since the reduction will depend on the form of the signals involved. However, from the 
discussion here it is apparent that the actual correlation function observed will be smaller in am- 
plitude than that predicted using an ideal image. This result is obvious if it is considered that the 
presence of the phase terms in the actual image produce a mismatch between the signal and refer- 
ence and therefore reduce the correlation. The effects of these phase terms can be reduced by 
further restricting the frequency range p max (or r max ), and signal, and image aperture size which 


45 



would limit the variation of the phase terms. An analysis of the required limitations will not be 
begun since it will depend to a large extent on the type of signals involved and the correlation re- 
sults desired. Here, the equations necessary for such an evaluation have been developed and hope- 
fully the significance of the various effects that appear in an optical system have been pointed out. 


PHASE CORRECTIONS 

In the last few sections were discussed the effects of undesirable phase terms in optical sys- 
tems and it was demonstrated that these effects can be minimized by restricting the size of signal 
apertures and the spectral range of the signals. An alternative approach can be pursued by insert- 
ing phase corrections into the optical system. Such phase corrections can be implemented by in- 
serting sheets or plates of transparent materials whose thickness or index of refraction has varia- 
tions that introduce phase terms opposite in sign to those introduced by the system. 

The basic equation representing the Fourier transform operation of a lens was given by Equa- 
tion 45 as 
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Rewriting this equation and retaining only the variable part of the terms outside the integral 

A(x, y) = <?!>(x, y) JJa' (x t , y t ) 6 12 (PXl+qyi) dx x dy x , 

where 

_ i 2 kf[(l + r2/2f 2 )(l + r2/f2) -1 / 2 ] 

<p(x, y) - e , assuming d - f , 

as derived in Equation 119. As pointed out in all discussions above, the phase term <£(x, y) de- 
stroys the simple Fourier transform representation of lens focussing properties since the integral 
part of Equation 131 corresponds to a Fourier transform by itself. Thus consider inserting a phase 
correction plate into the back focal plane of a lens with transmission properties given by 


P(x, y) - A 0 e iG <t>* (x, y) , 


(132) 


where 


A q - constant , 
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constant , 


and 


,* x _ -i2kf[< l+r 2 /2f 2 )(l+r2/f 2 )" 1/2 ] 

<P (X, y) - e 

The light amplitude distribution appearing at the output side of the plate will be 

A(x, y ) P(x, y) = JJa* (x,, yj e l 2 w ( px i +cly 1 ) dXj dy t , ( 133 ) 

where the constant term A 0 e ic has been dropped and &(x, y) 4 ? (x, y) = 1. Thus by inserting a 
phase plate with transmission properties given by Equation 132 in the back focal plane of each lens 
in an optical system, the phase terms are eliminated. From the definition given by Equation 132 
it is found that the phase correction depends on the focal length f of the lens and the wavelength 
\(k = 2n/k) of the light. The phase correction is not dependent on the signal used, and therefore 
a phase plate can be made for the lens and wavelength to be used in the system. Of course, the 
correction of phase by this method requires an accurate technique for producing the phase plate 
and positioning the plate in the optical system. In any case, it has been shown that the elimination 
of undesirable phase terms is possible at least in theory. Any inaccuracies in production or loca- 
tion of the phase plate may be acceptable as long as the phase terms are appreciably less than 
before the plate was introduced. Assuming the phase plate is an accurate representation of the 
transmission function of Equation 132, the Fourier transform relation of Equation 133 will be valid. 
With the relationship given by Equation 133, the operation of spectrum-analyzer, imaging, and 
optical- correlation systems can be described by the ideal cases used in the respective discussions, 
and no undesirable phase terms appear in the equations. 


PHASE TERMS WHEN d t f 

In the consideration of phase terms was considered the special case of an input plane coincident 
with the front focal plane of a lens (d = f). This special case was chosen to eliminate the phase 
effects of a term proportional to (f - d). From Equation 118 can be gotten a complete expression 
for the variable part of the exponential term e ikE as 


ik(d-f )( l+r 2 /f 2 ) 1/2 i2kf ( 1+ r 2 /2f 2 )( 1 + rVf 2 )“ 1/2 

<£(x, y) - e e 


(134) 


Equation 134 reduces to the form of Equation 119 when d = f . Restricting our consideration to a 
rather limited range in the back focal plane of a lens, it has been shown that Equation 134 can be 
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found to a good approximation in the form of Equation 88 such that 


4(x, y) = e ik(f-d)(r/f) 2/2 e ikf<r/f) 4 /4 _ (135) 

Equation 135 can' be considered as a representation of the phase in a back focal plane containing a 
frequency spectrum, while Equation 134 is a more accurate representation that applies in a back 
focal plane containing an image of the input signal. This application of Equation 134 and 135 is 
based on the relatively larger apertures commonly used in the signal and image planes. 

In the systems that have been considered here, the complete phase variations, as given by Equa- 
tion 134, appear only in the correlator system. This can be seen by noting the presence of <£(x 2 , y 2 ) 
in Equation 120. Since this phase factor is expressed in terms of the coordinates x 2 and y 2 of the 
image phase, a very restrictive aperture limitation cannot be used without severaly affecting our 
signal handling capability. Therefore, the approximation of Equation 135 will not be valid and 
<£(x 2 , y 2 ) in Equation 120 will have the form of Equation 134 with d 2 and r 2 replacing. d and r 
respectively: 
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where r 2 2 = x 2 2 +y 2 2 and d 2 is the distance from the spectrum plane F 1 to the lens L 2 (see Figure 
16). In the sample correlation function of Equation 130 it can be seen that the additional phase 
term dependent on (d 2 - f ) will increase the effect of <£(x 2 ) on the integrals. In practice, a system 
would be specified on the basis of locating lens L 2 so that d 2 = f. However, the exact positioning 
of the lenses in an optical system is obviously a practical impossibility. Thus the additional phase 
term containing (d - f ) represents the phase distortion introduced by inaccuracies in the imple- 
mentation of the system. Since the quantity ( d - f ) represents an inaccuracy, its value will usually 
be undetermined. Therefore, the first term in Equation 136 represents an undetermined phase 
error in the optical correlator system. If a guess or estimate of the tolerances in tl^e system can 
be made, this error term can be used to determine the maximum distortion of the correlation func- 
tion by analysis similar to that implied by Equation 130. 

Since the variation of the phase term containing (d 2 - f ) in Equation 136 is not known specifi- 
cally, the elimination of this term by a phase correction plate is not possible. Thus in a system 
containing phase correction plates, only the second term of Equation 136 can be eliminated. In 
such systems, 4 >{^ 2 , y 2 ) is completely given by the position error term 


^( x 2’ y 2 ) 
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(137) 


The distortion of the correlation function in a phase corrected system is, therefore, completely 
dependent upon the positioning errors. Again referring to the sample of Equation 130, <£(x 2 ) would 
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be given in the form of Equation 137. The phase term in front of the integral of Equation 130 would 
also be replaced by an error term from an expression such as Equation 135 as will be discussed 
below. 

The phase term given by Equation 135 represents the variable part of Equation 88. The ex- 
ponential dependent on (f -d x ) represents an error term due to lens positioning. To account for 
this error the complete phase approximation of Equation 135 must be used in place of the K : expo- 
nential of Equation 89. Thus the error phase term will appear throughout the previous analysis 
wherever the term has been used. 

The effect of the phase term has been considered on the image intensity and on the correla- 
tion function in earlier sections. In the correlator discussion, the variable phase term of K x appears 
in the integral used to define the image amplitude distribution in Equation 120. To account for 
errors in placement of lens Lj (see Figure 16) the exponential e iltf ( r / f ) 4 ^ in Equation 120 must 
be replaced by a phase term of the form of Equation 135, which can be written 
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Referring to the sample correlation of Equation 130, the phase term given by Equation 138 will 
replace the e ikf ( r - /f ) 4/ 4 term in front of the integral. The error phase term has the affect of adding 
an additional variation to the phase of the terms in the summation. For a phase corrected system, 
the e ikf c/ f > 4/ « term is eliminated and the undesirable phase difference of terms in the summation 
will be dependent only on the accuracy of the system implementation. 

In the discussion of imaging systems, a factor m was defined by Equation 103 and a method for 
determining the accuracy of the image intensity was developed based on this parameter. To extend 
this method to include the case for d / f , m can merely be redefined by 

m4 = 1 . (?) [( f_ d i) + 1 (?) ] > (139) 

which is obtained from the exponents of the terms in Equation 138. For d i less than f, the limits 
on m defined in the previous discussion will apply to Equation 139 for the maximum value of r. It 
is noted that since (f -dj is a positive quantity when d ( is less than f, the required limit on r max /f 
will be less than that determined for the case d ; = f . When d is greater than f , (f - d t ) is a nega- 
tive quantity which would imply that the value of r max /f can be greater than that for the case d = f . 
This is true except for cases in which d is sufficiently greater than f so that for some value of 
r/f less than r max /f the value of Equation 139 is greater in absolute value than for r max /f. That is, 
since (f -d t ) is negative, the right side of Equation 139 is zero at r = 0 and becomes negative as 
r increases until it reaches a maximum negative value and then increases to positive values. De- 
pending on the value of d t and the limit r nax /f , it is possible that the phase at the maximum nega- 
tive value is greater than that at the aperture limit r max /f. In such cases, the maximum negative 
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value must be considered rather than the end value at r max /f , and the aperture may have to be re- 
stricted to values below this maximum. In any case, the brackets on the right side of Equation 139 
must be considered as an absolute-value symbol when the quantity within is negative so that m will 
have real values. In other words, the concern is with the magnitude of the phase variations and 
not the sign. 

For systems with phase correction, only the first term of Equation 138 will remain and Equa- 
tion 139 will be simplified to 


m 


4 



(140) 


Equation 140 can then be used with the image intensity criterion developed earlier to consider the 
effects of positioning error for lens L x in phase corrected systems. 

In most of the literature the phase corrected form of Equation 138, 
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is used even though phase correction techniques may not be employed. This application of Equa- 
tion 141 requires that the frequency limitation be sufficient so that the e ikf c/ f ) 4/ 4 term can be 
neglected. This application also implies that the term (f -d x ) is much greater than the maximum 
value of f/2 ( r/f) 2 . If this condition does not hold, the neglected term will contribute a phase com- 
parable to that of Equation 141 which would then be in error. Conversely, if (f - d x ) is not greater 
than f/2 (r/f) 2 , and the e ikf < r / f > 4/ 4 term is considered negligible, than the term given by Equa- 
tion 141 is also negligible since it is comparable to the neglected term. 

In this section has been outlined the procedure for taking into account the additional phase 
term arising from inaccuracies in the positioning of lenses. It was pointed out that since these 
terms are caused by inaccuracies, they are generally not specified completely. The worse case, 
however, can be specified by estimating the maximum error in the position of a lens. From this 
extreme estimate the necessary aperture limitation, or the evaluation of errors in the desired op- 
tical outputs, can be determined for a worse-case analysis. Unfortunately, because of the unde- 
termined nature of these terms, phase correction cannot be used to eliminate their effects. 


SUMMARY 

The derivation presented in this report demonstrates that the Fourier transform representa- 
tion of a focussed diffraction pattern is a reasonable approximation for describing the operation 
of coherent optical systems with lenses. The basic assumptions consisted of the ideal focal proper- 
ties of a lens and the use of perfectly coherent light. Except for undesirable phase effects, it was 
demonstrated that the Fourier transform representation is obtained as a good approximation by 
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imposing limitations on the size and frequency content of signals allowed. The phase terms can 
also be eliminated by aperture limitations; however, the restrictions are more severe. Depending 
on the application, a trade-off must be made between the limitations required for elimination of 
undesired terms and the desired signal size and frequency content. 

Techniques for evaluating the effects of the various approximations and for analyzing the opera- 
tion of ideal optical systems have been presented. For specified signals and applications these 
expressions can be used to determine the theoretical errors in assuming ideal operations as is 
commonly practiced, The analysis presented is by no means complete; however, it is hoped that 
it is sufficiently detailed to provide a clear insight into the required approximations. 

This report represents an initial step in the development of a detailed analysis of the capabilities 
of optical processing systems. Further studies are required to formulate complete criteria and 
analysis techniques for practical optical systems. Some of the important areas which must be 
considered include: 

1. Lens aberrations 

2. Coherence 

3. Transmission properties of modulation media 

4. Band- limited signal approximations 

These areas were not treated in the analysis presented here since the initial study was restricted 
to ideal systems. The complexity of the mathematical formulation of optical patterns can be sim- 
plified somewhat by using the notation of communication theory (Reference 5). Such methods are 
becoming quite useful in modern optics studies. The development of these techniques provides a 
means for avoiding the complicated mathematical formulations inherent in diffraction problems. 
However, any new formulations such as these must be considered in terms of the more rigorous 
formulation since the various approximations are basically the same in both formulations. 


Goddard Space Flight Center 
National Aeronautics and S,pace Administration 
Greenbelt, Maryland, March 28, 1969 
125-23-02-05-51 
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Appendix A 


Derivation of Diffraction Formula 

This discussion is restricted to the somewhat special case of diffraction at an infinite plane 
surface as diagrammed in Figure Al. The shaded area in the figure represents a cross-section 
of an infinite slab. The basic problem is to determine the electric field at any point p in the dif- 
fraction region, which includes all points to the right of the plane boundary surface indicated in 
Figure Al. When the plane slab is not present, the electric field at any point P can be found simply 
by substituting the coordinates of P into the mathematical expression describing the light propagat- 
ing from whatever light source may be present. In itself, finding a mathematical representation 
for a given light source is not a simple problem. The light radiated by a square is dependent upon 
the mechanism generating the light as well as the geometry of the source. In many cases it is as- 
sumed that a good approximation is obtained by considering ideal light sources which radiate spher- 
ical waves (point sources) or plane waves (point sources at an infinite distance). 

Inserting a plane surface into the path of the light waves as shown in Figure Al complicates 
our problem. Since the presence of the plane effects the propagation characteristics in space, the 
electric field at any point will now depend on the characteristics of both the light source and the 
plane. The characteristics of the plane depend on the type of material of which it is made and these 
characteristics usually vary from point to point in the surface. Thus the problem is that of de- 
termining the electric field in the presence of a surface which can have widely varying electrical 
properties from point to point. As the reader 
may already know, problems of this nature are 
very difficult and, in fact, very few diffraction 
problems have been solved rigorously. For- 
tunately, in many cases of practical interest, 
results within experimental accuracy can be 
obtained by less rigorous techniques. 

In order to implement a discussion of dif- 
fraction problems, a formula must be derived 
for diffraction at a plane surface. This result 
was first derived by Sommerfeld 6 in 1896 and, 
as will be demonstrated, is effectively a math- 
ematical representation of Huygen's principle 
for the special case of a plane diffraction sur- 
face. The basic assumption we will start with 
is that the components of the electric field are 
known at every point on the right hand boundary 
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Figure Al— Diagram of diffraction configuration. 
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surface of the plane slab (refer to Figure Al). 
The methods for determining these field values 
are of no importance at this point; however, in 
many cases of interest, the assumption of a 
multiplying factor representing the transmission 
properties of a thin material provides results 
in close agreement with experiment. For the 
present purpose, it is simply assumed that the 
value of the electric field at every point on the 
plane boundary surface is known (i. e. , can be 
found easily). 


Referring to Figure A2, the problem to be 
solved can be stated as follows: 


Given the electric field at every point on an 
Figure A2— Outline of diffraction problem. infinite plane boundary, what is the electric 

field at any point p in the diffraction region? 


As shown in Figure A2, the diffraction region is defined to be all space on the right side of the 
boundary surface (note that this region does not contain any light sources). In Figure A2, it is as- 
sumed that all light sources are to the left of the boundary surface and that the diffraction region 
includes all points to the right of the boundary. It can be assumed that the diffraction region is in 
free space (velocity of light is c in all directions) and that there are no electric currents or charges 
present in this region. Since the electric field is a vector quantity, its direction at any point is as 
important as its magnitude. As in the case of any vector, the electric field can be considered in 
terms of its components in the x, y, and z directions. To simplify the discussion it is assumed that 
the light waves are monochromatic, or are constant in time at a single frequency. When necessary 
this discussion can be extended to the general case of nonmonochromatic waves by considering each 
separate frequency component, as described here, and summing up all components. 


Each component E (in x, y, z components) of the electric field of a monochromatic wave will 
satisfy the Helmhotz equation (time -independent wave equation) at every point p in free space which 
contains no electrical sources: 


( V 2 +k 2 )E = 0 , 


(Al) 


where 


E - x, y, or z component of the electric field ( 


d 2 d 2 + d 2 
x 2 dy 2 <9z 2 
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CO 2 TT 

k - c = X ’ 

co - angular frequency , 
c - speed of light , 

K - wavelength of light 


Using the values for e on the plane boundary surface and the fact that e must satisfy Equation A1 
at every point P in our diffraction region, a formula can be derived for the electric field at P in 
terms of the values given on the surface. 

An arbitrary function V which also satisfies Helmholtz equation can be introduced 

(v 2 + k 2 )v = 0 . (A2) 

There are many functions which will satisfy Equation A2; however, we will continue to use the sym- 
bol v and reserve the selection of a specific function until a few additional characteristics of v 
which will allow us to accomplish our derivation have been determined. In terms of E and v, two 
vectors F and f 2 can be defined as 


and 


Fj = E V V , 


(A3) 


F 2 = V VE , 


(A4) 


where 


V 


<5 . v d - 5 
1 + 3 + k al • 


and W and VE denote the gradient of v and e respectively. Gauss's theorem is introduced as 


jJJv-Fdv 

volume 


|J ? ■ di 

su r face 


(A5) 


where the volume integral on the left can be taken over any volume that does not contain discon- 
tinuities of the divergence of F(V • F) and the surface integral on the right is over the surface which 
encloses the volume (f must be continuous on the surface so that the surface integral can be found). 
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From Equations A3 and A4, 


V • F = V • E V V = VE • W + E V 2 V 


(A6) 


and 


V • F 2 = V • WE = VV • VE + VV 2 E . 

From above it is noted that Equations A6 and A7 must not have discontinuities within the volume of 
integration; therefore, E and v must have continuous first and second derivatives. Since diffraction 
region free of electrical current and charge is under consideration, E will meet this requirement 
for any volume in the diffraction region. Since a specific v has not yet been selected, we will note 
this requirement and be sure to satisfy it when selecting v. Thus Equation 5A can be written by 
substituting Equations 3A, 4A, 6A, and 7A: 


III 7 • *■ * = JJJ{VE • W + EV 2 v}dv = {Jew . ds 


JJJv • f 2 dv = JJ{{vv • VE+ W 2 E}dv = {{ 


VVE • ds 


(A8) 


(A 9) 


Equation A9 can be subtracted from Equation A8, and noting that VE • W - W • Ve = 0, then 


{{{{eV 2 V- W 2 e} dv = {{ 


{EW-WE} • ds 


(A10) 


This is Green’s theorem and holds for any functions E and v which have continuous first and second 
derivations in the region of integration. From Equations A1 and A2 it is known that v 2 e = -k 2 e 
and V 2 v = - k 2 v and therefore the bracket on the left hand side of Equation A10 gives 

EV 2 V-W 2 E = - k 2 EV - (-k 2 EV) = 0 . (All) 

Since the integrand is zero as given in Equation All, the volume integral on the left hand side of 
Equation A10 is zero and Equation A10 can be rewritten as 


I 1 


{ EW - WE > • ds 


(A12) 
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The surface integral in Equation A12 is to 
be taken over any closed surface which does 
not enclose discontinuous points. For the pur- 
poses here, the surface is chosen as indicated 
in Figure A3. The outer surface consists of the 
infinite boundary plane, A, on which the electric 
field is known, and a hemisphere, C, of infinite 
radius which connects the ends of the plane at 
infinity. An inner surface,!, is defined as a 
sphere centered at P with radius e . Taking A, 
C, and ! as the closed surface, the volume be- 
tween the sphere ! and the outer surface, A - C, 
is defined. If the limit is taken as e goes to 
zero, the point P will be the only point in the 
diffraction region outside the surface of inte- 
gration. Thus the sphere ! isolates the point P 
where we want to find the electric field. The 
surface integral of Equation A12 can be written 
as the sum of the integrals over A, c, and !: 


PLANE BOUNDARY 

SURFACE, A 

(ELECTRIC FIELD KNOWN 


CO 

Figure A3— Boundary surface enclosing all points except P. 



II 


JJ< 


{EW - WE} * ds + {EW - WE} * ds + 


* 


{EW - VVE} • ds 


(A13) 


The surface integral over the hemisphere at infinity can be eliminated through a physical ar- 
gument given by Born and Wolf (Reference 4). In practice a light wave starts at some time, and 
since it propagates at a finite velocity (c in free space) must have an end. We can imagine the in- 
finite hemisphere continually expanding in front of the light waves. In this way, the contribution of 
the wave on the hemisphere is zero since the light waves never reach the hemisphere. The integral 
over the surface c will therefore be zero and Equation A13 can be written 


II 


(EW - WE} ‘ ds + {EW - WE} • ds 


II 


(A14) 


Advantage is now taken of our freedom to select a function v in order to simplify Equation A14. In 
the integral over the plane surface A it is noted that the term WE requires the values of VE on the 
boundary surface. Since VE is not known on A, v must be required to be zero on the boundary surface 
to eliminate this term. Equation A14 can then be written 


JJeW ■ ds + JJ- 


{EW-WE} • ds 


(A 15) 
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To consider the integral over the sphere £, the surface element ds is expressed in polar coordinates 

ds - “ € 2 s in # d# dc£ r , 


where f is a unit vector in a radial direction away from the center at p and the minus sign is re- 
quired by the convention that a surface normal is directed away from the enclosed volume. The 
integral over £ can be written as 


- r. ta »* r 

Jo Jo 


d<£ € 2 (EW~ WE) * ? , 


(A16) 


By vector identity (EW- WE) • f = E (dV/d r) - v (dE/dr), and Equation A16 can be written 


r 

r r 

/ , dv) 

1 o 

1 s in 9 dO 

h 

l \ e dr J 

-fJ 


Jo 




(A17) 


But the desire is to take the limit as e goes to zero. The term e 2 v (dE/dr) in Equation A17 can be 
eliminated if v is required to satisfy the condition 


lim e 2 V - 0 , 

e-0 


(A 18) 


at any point P in the diffraction region. 

Since e has continuous first derivations <9E/dr will be finite and the conditioning given by Equa- 
tion A18 will give 


Expression A17 is then given as 


lim 

e-o 


dE 

v d7 = 


(A19) 


£ sin0d0 E ( fi2 37 )} 


(A20) 


Since E is the value of the field on the surface £, and the surface £ reduces to the point p when e 
goes to zero, the limit in Equation A20 can be written 


lim Ee 2 
e-0 


dV 

dr 


E(P) lim e 2 
£-0 


dV 

dr 


(A21) 
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where E(P) is the electric field atP. We will require that V satisfy the condition 


e2 TE = 1 . (A22) 

at any point P in the diffraction region. The limit in the brackets of Equation A20 is then simply 
E(P) and Equation A20 gives 


- r,i„ 9 *r 

•'0 Jo 


sinSdS d£E(P) = - 4 ttE(P) . 

'o 


(A23) 


Substituting the result of Equation A23 into Equation A15 for the surface integral over Z, 


S 


EVV • ds - 4 ttE(P) = 0 


(A24) 


Rearranging terms, Equation A24 can be written 


E(P) 


= i JJ 


EW • ds . 


(A25) 


By vector identity, W • ds = (dV/dn) ds , where 3 / 3 n represents a partial derivative with respect to 
a coordinate perpendicular to the plane boundary surface in a direction out from the enclosed re- 
gion. Equation A25 can be written without vectors as 


E(P) 



(A26) 


Except for the selection of a function v which satisfies all the conditions that have been used here, 
Equation A26 has the required form. The left side is just the field at a point P, and since the inte- 
gral on the right side is on the plane boundary surface, a the e in the integrand assume the given 
values on the surface. Thus the field at any point p is given in terms of the given values on the 
surface A by Equation A26. 

In deriving Equation A26, restrictions have been imposed which must be satisfied by the func- 
tion v. The requirements to be satisfied by v are: 

1. V must have continuous first and second derivations within the region inside the boundaries 
shown in Figure A3. 

2. (V 1 2 + k 2 ) V = 0. 
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3. V = 0 on boundary surface. 

4. W £ 0 on boundary surface. 

5. iim e 2 v = 0 at any point p in diffraction region. 

e -»0 

6. Iim e 2 (dV/d r) = 1 at any point p in diffraction region. 

e-0 

Fortunately there is a function which meets all these requirements: 



(A27) 



Figure A4— Geometry for definition of V. 


where r and r ' are defined by Figure A4 
as 

r = distance from P to any point Q , 
r ' ~ distance from P' to any point Q, 


P' 


mirror image of P 

(i.e., PP' is perpendicular to 

the boundary A and d - d' ) , 


and 


co _ 2n 

c ~ k 


That the function v given by Equa- 
tion A22 does in fact satisfy all require- 
ments, can be proved by direct substi- 
tution into the expressions listed above. Here we will only discuss the continuity requirements 
(item 1 above). As defined by Equation A27, V has discontinuities at r = 0 and at r ' = 0. These 
discontinuities appear at the points p and P' respectively; p' lies outside the diffraction region, 
andp was separated out of the integration region by the sphere 2. Thus the only points at which 
discontinuities appear are outside the region specified in the continuity requirement, and v given 
by Equation A27 does satisfy this requirement. 


Returning to Equation A26, it is noted that dV/dn on the surface A is required, rather than V 
itself. Coordinates can be selected so that the z axis is perpendicular to surface a as shown in 
Figure A5 so that dV/dn becomes - (dV/dz) (minus sign appears since the positive Z direction is 
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opposite to the positive n direction). From 
the geometry of Figure A5, r and r ' are given 
as 

r = [x 2 +y 2 + (d “ z) 2 ] 1/2 , (A28) 

and 

r* = [ X 2 +y 2 +(d + z )2]V2 . (A29) 

Substituting Equations A28 and A29 into Equa- 
tion A27 ? it is found that 



Figure A5— -Geometry for definition of terms in V. 


(9V <9V S fe ik ^ x2+ y 2+ ( d+z ) 2 ] 1/2 e ik[x 2 +y 2 +(d-z) 2 ] 1/2> l 

^ Z ifx 2 + Y 2 + (d+ z) 2 j 1/2 [x 2 +y 2 + (d - z) 2 ] 1/2 j 

= |---- ik ( d+ Z ) + 1 e ik[x 2 +y 2 +(d+z) 2] 1/2 

|Jx 2 +y 2 +(d+ z) 2 ] [x 2 +y 2 +(d+ z) 2 ] 3/2 J 


-ik(d~ z) 


[x 2 +y 2 + (d - z) 2 J [x 2 +y 2 + (d- z) 2 ] 


3/2 


i k Cx 2 +y 2 +(d~z ) 


2 ] 1/2 


(A30) 


Now in Equation A26, integration was over the surface a, so that the value of - (dV/dz) on A which 
is obtained by setting z = 0 in Equation A30 is needed: 


dV 


A 


dV 

dz 


z-0 


_ 2 e iktc 2 + y 2 + d 2 ] 1/2 d 

[x 2 +y 2 + d 2 ] 1/ 2 [ x 2 + y 2 + d 2j 1/2 



1_ 

[x 2 + y 2 + 



(A31) 


Equation A31 can be simplified by noting from Equation A28 and the geometry of Figure A5* that 


r - [x 2 + y 2 + d 2 ] 1/2 , when z-0 


(A3 2) 


cos 0 


d 

[x 2 + y 2 + d 2 ] 1/2 


(A3 3) 


*Note that in Figure A5, P and P # were chosen as point on the z axis to obtain Equations A32 and A33* In general x would be replaced 
by (x - x Q ) and y would be replaced by (y - y Q ) where x Q , y Q define the x, y coordinates of the points P and P \ In the derivation, 
the z coordinate of the point P was represented by d to avoid confusion with the coordinates of the point Q. In general the d in Equa- 
tions A32 and A33 is replaced by z. Otherwise the general results have the same form as found above. 
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Substituting Equations A3 2 and A3 3 into Equation A31, then 


dv 

dn 


A 


dV\ 

dz\ 


2 = 0 



(A3 4) 


Now Equation A34 can be substituted into Equation A26 so that the diffraction formula can be com- 
pleted as 


E(P) 



(A3 5) 


where (refer to Figure A5) 

E(P) = Electric field at a point p in the diffraction region, 
e a = electric field on the plane boundary A, 
r = distance from P to a point on A, 

6 = angle between r and normal to plane A, 


and 


, _ co 2v 
k = c = X • 

Given the electric field e at every point on the surface A, Equation A3 5 can be used to determine 
the electric field E(P) at any point P in the diffraction region. This statement depends, of course, 
on whether the integration indicated in Equation A3 5 can be performed. If the integration cannot 
be performed analytically, it can be assumed that a numerical solution to any desired accuracy can 
be obtained using a computer. In many problems of interest, satisfactory results can be obtained 
by approximating Equation A35 using the geometry of the specific problem. For example, for small 
angles d such that cos 9 ~ 1, and at great distances r such that 1/r « k, Equation A35 can be ap- 
proximated as 


E(P ) 



(A3 6) 


Equation A36 represents Huygens' principle since the contribution from each point on the boundary 
surface is given by e a e ikr /r which describes the spatial variation of a spherical wave. Thus the 
diffracted field as given by Equation A36 can be interpreted as the summation (integral) of spherical 
waves radiating from each point on the diffraction boundary. 
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